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ABSTRACT 


The Electrical and Computer Engineering Department at the Naval Postgrad- 
uate School has a need for additional software to be used in instructing students 
studying signal processing. This software will be used in a PC lab or at home. 

This thesis provides a set of disks written in APL (A Programming Language) 
which allows the user to input arbitrary signals from a disk, to perform various 
signal processing operations, to plot the results, and to save them without the need 
for complicated programming. 

The software is in the form of a digital signal processing “toolkit.” The user 
can select functions which can operate on the signals and interactively apply them 
in any order. The user can also easily develop new functions and include them in 
the “toolkit.” 

The thesis includes brief discussions about the library workspaces, a user man- 
ual, function listings with examples of their use, and an application paper. The 


software is modular and can be expanded by adding additional sets of functions. 
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I. INTRODUCTION 


The Electrical and Computer Engineering Department at the Naval Postgrad- 
uate School has a need for signal processing software packages for use by students 
studying digital signal processing (DSP). The purpose of this thesis is to provide a 
set of tools in the form of APL workspaces on floppy disks that students can use at 
home or at school on IBM-PC-compatible computers for solving digital signal pro- 
cessing problems. Workspace is a collection of functions and data that are placed 
together in order to do any job. 

The software supports the following currently offered courses in DSP: 

1. EC 2400-Discrete Systems 

2. EC 3400-Introduction to Digital Signal Processing 

3. EC 3410-Introduction to Discrete Time Random reus 

4, EC 4430-Advanced Signal Processing and Spectral Estimation 
5. EC 4440-Multidimensional Digital Signal Processing 

Specific APL workspaces were developed for EC3400, EC3410, and EC4440. 
The software uses STSC Corporation’s APL*PLUS PC and STATGRAPHICS sys- 
tems which are both site licensed to the Naval Postgraduate School. APL*PLUS 
PC is STSC’s version of the APL language for the IBM-PC and compatibles. Sim- 
ilar versions of the APL language are also available from STSC for other com- 
puters. STATGRAPHICS (Statistical Graphics System) is a PC software package 
integrating a variety of statistical functions with color graphics. STATGRAPH- 
ICS is written in APL (specifically APL*PLUS PC) and offers direct access to the 


host language. Currently STATGRAPHICS is available only for the IBM-PC and 
compatibles. 

The DSP software developed in this thesis is in the form of a digital signal 
processing “toolkit.” The user can select functions to operate on the signals and 
interactively apply them in any order. The user can also easily develop new func- 
tions and include them in the “toolkit.” The tools are provided on a set of floppy 


disks that users can carry to their own PC machine. The sets of disks contain: 
1. The APL*PLUS Executive Program 
2. STATGRAPHICS disks 
3. UTILITY workspace 
4. EC3400 workspace 
5. EC3410 workspace 


6. EC4440 workspace 
The main products of the research are the disks that include the signal pro- 
cessing library workspace, and the documentation provided in this thesis. A typical 
library contains functions such as: Fast Fourier Transform (FFT), Inverse FFT, 
Sine, Cosine, Sample Autocorrelation, Mean, and other similar functions. In de- 


veloping the DSP software the following requirements were met: 


1. Ability to input/output from/to a disk data file in free format. This allows 
convenient communication with DOS and with programs written in other lan- 
guages or on other computer systems. 


2. Ability to provide for graphic output to represent signals, frequency responses, 
and other functions including 2-D, 3-D, and contour plotting. 


3. Ability to perform filtering, convolution, and other common signal processing 
operations in the signal (time/space) domain. 


4. Ability to compute transfer functions, to compute frequency responses, and 
generally to perform complex arithmetic and other common frequency domain 
operations. 


5. Ability to generate random signals. 


6. Ability to compute statistical characteristics of the signals including correla- 
tion functions, spectra, and prediction errors. 


The user can combine functions in any desired order to solve DSP computer 
assignments and can concentrate on the signal processing exercises without the 
need to do complicated programming. 

Although the intent was to develop the signal processing toolkit specifically for 
the IBM-PC and compatibles, the system is portable. STSC has APL available on 
IBM mainframes, Digital Equipment Corporation VAX, (under UNIX and VMS), 
Apple MacIntosh, and others. These DSP workspaces are portable to all of these 
systems with few, if any, changes. In addition, the DSP workspaces are easily 
carried over to other APL systems with minor changes to only some of the most 
basic functions; bits and pieces of the software have existed for some time under 
VS-APL on the IBM mainframe and under UNIX 4.3 bsd. 

The remainder of the thesis is organized as follows. 

Chapter II contains a brief description of the APL language and the STAT- 
GRAPHICS software system. Chapter III deals with the DSP library and discusses 
the application of the four workspaces. The Utility workspace includes general 
functions such as GETDATA and PUTDATA for the data input/output to DOS; 
the EC3400 workspace includes basic one-dimensional DSP functions; the EC3410 
workspace includes one-dimensional DSP functions for statistical signal processing; 
and the EC4440 workspace includes multidimensional DSP functions. Chapter IV 


outlines the conclusions, and the benefits of using the entire package as a software 


toolkit. Appendix A provides a simple user manual for the software package. Ap- 
pendix B provides a function library with a brief description of each function. Ap- 
pendix C contains function listings and samples of their use. Appendix D furnishes 
samples of computer assignments and their solutions using the tools developed in 


this thesis. 


II. APL AND STATGRAPHICS ENVIRONMENT 


The signal processing software is written in APL and uses the APL*PLUS PC 
interpreter and the STATGRAPHICS software system. Although STATGRAPH- 
ICS provides many statistical functions in itself, only the plotting functions from 
STATGRAPHICS are used. Following is a brief description of APL and the STAT- 
GRAPHICS system. 


A. APL [1] 

Originally APL began as a notation to express mathematical procedures. Its 
originator, Dr. Kenneth Iverson, published his notation in 1962 as A Program- 
ming Language from which the name “APL” is derived [2]. 

One of the strongest advantages of APL is its ability to handle complicated 
array data-structures without extensive programming. The language allows the 
definition of numerical items or a series of items with one symbol. Examples of 
data that can be represented are: 

e SCALAR - A single value 

e VECTOR - A group of scalars 
e MATRIX - A group of vectors 
e ARRAY - A group of matrices 

The language makes it as easy to perform an arithmetic operation using a 
matrix or vector argument as to perform the operation using simple numbers. 


This characteristic is illustrated in the following example. 


Μει6 (define vectors with 6 numbers) 


Ἱ ο το 


Ne2xv6 
N 
20496579 10 12 


NxM (vector multiplication) 
2 8 18 32 50 72 


Mes 3.799 
M 


Moo 
نہ ایا ی0‎ 
ιο σι ω 


(define matrices with 3 rows and 3 columns) 


Ne3 3 p3x19 
N 


3 6 9 
12 15 18 
20021157 
NxM 
351207 -- 
48 75 108 (matrix arithmetic operation) 
147 192 243 


Notice that when two items are of the same dimensions, the arithmetic operation is 
performed element by element. Thus N x M is not the usual linear algebra product 
of matrices. The latter is realized by a composite operation denoted by +. x. 


N+. xM 

90 108 126 
198 243 288 
306 378 450 


Scalar arguments are automatically extended to be conformable with the arrays. 


Thus the statement M x 2 yields 


2 4 6 
8 10 12 
14 16 18 


There is frequently no need for loops or similar mechanisms as in other lan- 
guages such as Fortran. From the user’s point of view the language performs the 
operation on the elements of the array in parallel without the need for loops. 

Note also that there is no need for dimensions because APL carries the dimen- 
sions of an item along with its definitions. 

APL provides a function to help the user to define or to determine the number 
of elements in the variable. The function represented by the Greek character “p”, 
is called “shape,” since 1t gives the form of the variable. On the previous page 
there is an example where we define the matrices M and N using the function “p”. 
In order to find the shape of any variable, use the function “p” with the variable 


name following it. The following is an example. 


pM 
3 3 


APL has a large variety of mathematical functions that can easily manipulate 
the data and can handle operations on a simple as well as on a complicated level. 
Thus, it is possible to place on a single line a combination of several functions which 
work together in a pipeline. All of this can be done interactively. The flexibility in 
manipulating array data easily makes the language well suited for single processing 


analysis. The following are a few examples: 


, جا‎ i EET (outer product) used to define separable 
9 2.0 7 2-D function z(n1,n2) = T(n lz 20a) 
10 MIS 
EN x3 
Me3 3 9135790864 
M 
mas 
790 
364 
The domino operator (calculate M inverse) 
gM 


"9.138313131313 70۰09090909091 0.2272727273 

Q. 1414141414 0. 1818181815 “0.1767676768 

0ϱ.1515151515  70.09090909091 0.06060606061 
Cells 6 


(solve Mz = C for z) 


CIM 
70.9090909091 1.04040404 1.757575758 


(check solution by 
multiplying constant vector 
by M inverse) 

T. X is matrix product 


(EM) +. XC 
70.9090909091 1.04040404 1.757575758 


APL uses symbols to represent its own built-in functions and, in addition, APL 
allows the user to define other functions. The solution to one problem consists of 
a collection of the programmer’s functions, each one of which can be applied inter- 
actively by the user or can be organized to call another function. The advantage 
of this organization is the ease of debugging and that only a few lines of APL can 
do volumes of work. 

APL has a concept called the workspace which is a collection of functions and 
data that are placed together in order to do a job. That collection is stored in the 


computer memory during an APL session and gives the operator immediate 


feedback and accessibility, and the workspace can be saved on disk and loaded back 
into the memory at any time for future use. 
B. STATGRAPHICS SOFTWARE PACKAGE [3] 

STATGRAPHICS is a Statistical Graphics System written in APL for scientific 
and engineering applications and for data analysis. The software supports dot 
matrix printers, monochrome/color terminals, and pen plotters. It runs under the 
APL*PLUS PC interpreter version 5.0 or later, which allows the user to write 
his own functions. Using STATGRAPHICS in this mode requires familiarity with 
APL. 

The system also has a menu that contains graphic functions, a set of data 
analysis and statistical procedures. This menu is illustrated in Figure 2.1. In 
addition, the menu contains data management and system utilities that allow the 
user to select the system environment, control color, plot dimensions, and so on. 
The most often used screens for DSP applications are: 

1. Plotting functions (Figure 2.2) 
2. Time Series Analysis (Figure 2.3) 
3. Experimental Design (Figure 2.4) 
STATGRAPHICS has options for modifying graphs and for on-line help. 

To use STATGRAPHICS with APL the user needs to start APL*PLUS PC as 
usual and then load the STATGRAPHICS workspace from the start-up disk. By 
going back and forth from the APL environment to STATGRAPHICS, the user 
can process data and plot the results easily. This makes the combination of APL 


and STATGRAPHICS powerful and very convenient to use. 


STATGRAPHICS Statistical Graphics System! 


DATA MANAGEMENT AND SYSTEM UTILITIES DINE SERIES PROCEDURES 
A. Data Management L. Forecasting 
B. System Environment M. Quality Control 


C. Report Writer and Graphics Replay N. Smootning 
D. Plotter Interface O. Time Series Analysis 


PLOTTING AND DESCRIPTIVE STATISTICS ADVANCED PROCEDURES 
E. Plotting Functions P. Categorical Data Analysis 
F. Descriptive Methods Q. Multivariate Metnods 
G. Estimation and Testing R. Nonparametric Methods 
H. Distribution Functions S. Sampling 
I. Exploratory Data Analysis T. Experimental Design 
ANOVA AND REGRESSION ANALYSIS MATHEMATICAL AND USZR PROCEDURES 


U. Mathematical Functions 


J. Analysis of Variance 
V. Supplementary Operations 


K. Regression Analysis 


Use cursor keys to highlight desired section. Then press ENTER. 
lHelp s Edi rE 3Savser 4Prtser 5 6Go 7Vars 8Cad 9Review 10Quit 


INPUT 11/18/87 13:08 STATGRAPHICS Vers. 2.1 


Figure 2.1. STATGRAPHICS Main Menu. 





l. X-Y Line and Scatterplots 
2. Multiple X-Y Plots 

3. X-Y-Z Line and Scattarplots 
4. Multiple X-Y-Z Plots 

5. Barcharts 

6. Piechares 

7. Component Line Charts 


Use cursor keys to highlight desired procedure. Then press ENTER. 
lHelp ZEGE 3Savscr 4Prascr 5 6Go 7Vars 8Cad 9Review lO0Quit 


INPUT 11/18/87 13:10 STATGRAPHICS Vers. 2.1 


Figure 2.2. The Plotting Function Menu. 
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| TIME SERIES ANALYSIS 


l. Horizontal Time Sequenca Plot 

2. Vertical Time Sequence Plot 

3. Seasonal Subseries Plot 

4. Autocorrelation Function 

5. Partial Autocorrelation Function 
6. Cross-Correlation Function 

7. Simple or Seasonal Differencing 
8. Mean or Trend Removal 

9. Box-Cox Transformation 

10. Periodogram 
11. Integrated Periodogram 

12. Data Tapering 
13.  Plotting vs. Fourier Frequencies 
14.  Box-Jenkins ARIMA Modeling 
15.  Cross-Correlation Matrix Plot 


Use cursor keys to highlight desired procedure. 


. 0 2598 3Savser 4Prtscr 5 6σο 7Vars 


INPUT 11/18/87 13:13 STATGRAPHICS Vers. 2.1 


Figure 2.3. The Time Series Analysis Menu. 


Then press ENTER. 
8Cnd 9Review 10Quit 


EXPERIMENTAL DESIGN 


1. Full and Fractional Factorials 


2. Central Composite Designs 
3. Alias Structure 
4. Response Surface Plotting 


Use cursor keys to highlight desired procedure. 
lHelo ταις 3Savscr 4Prtscr 5 6σο 7Vars 


Figure 2.4. The Experimental Design Menu. 
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Then press ENTER. 
8Cmd 9Review 1nOuit 


III. THE DSP WORKSPACES AND APPLICATIONS 


The following is the principle part of the thesis research. This chapter describes 


the workspaces that were developed by the author, and the way they can be applied. 


A. UTILITY WORKSPACE 

This workspace contains functions of general use that can be applied in combi- 
nation with functions from other workspaces in this package. Appendix B contains 
the list of workspace functions and a brief discussion of each one. 

Two of the most useful functions in the utilities are described and demon- 
strated below. The functions GETDATA and PUTDATA allow the user to transfer 
data from the APL environment to the DOS environment and vice versa. These 
functions use file system commands that are specific to STSC’s implementation on 
the IBM- PC. As such, the exact APL code for GETDATA and PUTDATA is not 
always portable to APL on other systems. However, versions of GETDATA and 
PUTDATA can be written for any system. Versions are currently available for the 
IBM mainframe under VS/APL and for a version of APL under UNIX 4.3 bsd on 
the VAX. To the user, GETDATA and PUTDATA provide a common interface for 
communication with other programs. 

The following is an example showing the transfer of data using the GET- 
DATA and PUTDATA functions. Figure 3.1 describes two sets of data in an APL 
workspace. In order to transfer those sets to a file in the DOS environment the 


user needs to type the following commands: 
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X PUTDATA 'X.DAT? 
Y PUTDATA ‘Y.DAT’ 
X is Š5 x 5 and Y is a 2 x 4 x 4 array 


x 
0.925 1.35 2.175! 3.7 4.625 
5.55 6.475 7.4 8.325 9.25 
10.175 11.1 1290251125990 13.379 
14.3 15.125.16.65 17.575 18.5 
19.425 20.35 21.275 22.2 23.125 


Y 
1£50 2E60 3560 4E60 
5560 6260 70 8E60 
9E60 2 “1544 “2E44 
W3E44 دو‎ “SE44 “6544 


7E44 -3E44 "9E44 “LE4S 

“1.1545 “L.2E45 “1.3845 "1.4E45 
WL.SE45 “1.6E45 “1.7E45 “1.3E45 
“1.9845 “2545 “2.1845 “2.2845 


Figure 3.1. X and Y Data Sets. 


The left argument to PUTDATA is the variable name; the right argument is the 
file name or path. Note that after کت‎ oe any data set to DOS the data is not 
erased and remains available in the workspace. 
When terminating APL, the user is placed in the DOS environment. Figure 
3.2 shows that the data set transfer was successfully completed. Note that the first 
line of the data set gives the dimension or “shape” of the data. 
By going back to APL environment and typing the M 
X—GETDATA 'X.DAT' 
Y—GETDATA 'Y.DAT? 
the data is transferred correctly from the DOS environment to the APL workspace. 
Figure 3.3 shows that the data has the same original shape. The data are read 
into the workspace and assigned to the variables X and Y. The data still remain 


available as a file in the DOS environment. 
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(terminate APL) 
OFF 


C:\APLoTYPE X.DAT (type of contents of files) 


5 5 
0.925 1.85 Le τσ QUY 
4.625 5.55 6.475 7.4 
8.325" 39225 FF 10. 1 ٦٣۶٠ 
12.025 12-35 89 
15.725 16:65 ۷٦ 
19.425 20.55, 1775۶ 
23.125 
C:\APLOTYPE Y.DAT 
2 4 4 
1E60 2E60 3E60 4E60 
5E60 6E60 7E60 8E60 
9E60 PEGI -l1E44 -2E44 
=3 E44 -4E44 -5E44 -6E44 
-7E44 -8544 -9E44 -LE45 


-1.1E45 -1.2E45 -1.3E45 -1.4E45 
-1.5E45 -1.6E45 -1.7E45 -1.8E45 
-1.9E45 -2E45 -2.1E45 -2.2E45 
C:\APL> 


Figure 3.2. And Y Data Sets In DOS Environment. 


ex 
5979 
X 
07925. 1.835 272759 7 4.625 
5.55 6.475 7.4 8.325 49.25 


575 1259901131 12.025 11.1 10.175 
5 76ء 16 "7257 ۔ 25 14.8 
1235 -23 22.2 20.3521200 19.42 


pY 
24 4 

T 
20 2E60 3E60 4E60 
5E60 6E60 7E60 8E60 
9560 LES “1E44 -2544 
“3544 "4E44 75E44 “6E44 
“7Ε44 “8E44 “9E44 “1E45 


"1 1845 “1.2E45 “1.3E45 "1.4E45 
71.5E45 -1.6E45 “21.7845 71.8E45 
"1.9E45 5 2. 1845. 2 E 45 


Figure 3.3. Getting X And Y Into APL Workspace. 


B. EC3400 WORKSPACE 
EC 3400 is a course at the Naval Postgraduate School that introduces first 


principles of digital signal processing. The topics covered include: the Discrete 
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Fourier Transform and the FFT algorithm, flow-graph and matrix representation 
of filters, ideal filters and approximation, and design of recursive and nonrecursive 
digital filters, fast convolution and correlation. 

The following workspace enables the student to solve computer assignments on 
the above subjects in a simple way. This workspace contains functions such as: lin- 
ear convolution, circular convolution, FFT, inverse FFT, and complex arithmetic. 
Some of those functions will be described later. In addition, the EC3400 workspace 
contains functions which allows the student to design one-dimensional filters. The 
methods used are the Fourier method and the frequency sampling method and 
windowing of the filters coefficients. The full description of the workspace is given 


in Appendix B. The following is an example of how to use this workspace. 


EC 3400 Computer Assignment {4 


| 
A highpass filter is to have an analog cutoff frequency of 6kHz. The filter is 


¡to be implemented by a digital filter having a sampling frequency of 40kHz. Plot x 


the frequency response of the filter with and without a Hamming window if 51 


coefficients are used. 





Solution 
The solution given here uses functions provided in the EC3400 workspace. In 
an actual classroom assignment, students may be asked to generate their own func- 
tions or to solve the problem in their own way. By using the function DIGFREQ, 
the digital cutoff frequency is calculated. (Figure 3.4). This function, like several 
others is included for convenience only; the calculation in converting radian to dig- 
ital frequency can be done easily in the workspace by using the built-in arithmetic 


operations. Figure 3.5 describes the ideal frequency response to be generated. 
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By using the function HPCOEFF, 51 causal impulse response coefficients are cal- 
culated. (See Figures 3.6 and 3.7). The function HPCOEFF uses the following 


equation: 


hren) = EAS al ap oH 


TC#40000 DIGFREQ 8000 
TC 
1.24665 7061 


TC+PI 
0.4 


Figure 3.4. Calculating a Digital Cutoff Frequency. 


H(e) 





ao 


theta ExPIr64I] 


Figure 3.5. The Ideal Frequency Response. 
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ΕΙ ECG ErS TC 


hHP 

6.8632316738E^17 70.01261377881 78.134689424E^3 8.504448034E^73 0.01441574 
2.239861191E 17 70.01593319429 ”0.01039432538 0.01100575628 
0.01892066822 2.2339861191E^17 ”0.02162362082 ”0.01439214283 
0.01559148806 0.02752097195 2.339861191E^17 70.0336367435 “ο. 02138 
0.02672826525 0.05045511524 2.339861191E^17 70.07568267286 
~0.06236595225 0.09354892838 0.3027306915 “0.6 0.3027306915 
0.09354892838 70.06236595225 70.07568267286 2.339861191£”17 
0.05045511524 0.02672826525 70.02338723209 “0.0336367435 2.3398611 
0.02752097195 0.01559148806 70.01439214283 ^70.02162362082 
2.339861191E^17 0.01892066822 0.01100575628 70.01039432538 
70.01593319429 2.339861191E^717 0.01441574721 8.504448034E73 
~8.134689424E°3 “0.01261377881 6.363316738E717 


Figure 3.6. Calculating the Impulse Response Coefficients. 


Sl coeff. for causal rnonrecursi ve 
HP filter 
O. 4 
O. 2 SE 
hin) o 


EC. = 





samples Em 
Figure 3.7. The Impulse Response Coefficients. 
The function HAMMING generates “J” samples of a shifted Hamming window 


using the following equation: 


2 
W(n) = 0.54 + 0.46 cos (735) πρ 


This is multiplied by the filter coefficients to obtain the windowed coefficients. (See 


Figures 3.8 and 3.9). 
By using the function FREQRES, the frequency response with and without a 


window is generated. This is done by calculating the DFT of the filter coefficients 
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WeHAMMING 51 


W 


0.08087246878 0.08783237415 0.1016466798 0.122105975 0.1489001176 0.1816 
0.31 0.3606984983 0.4141150246 0.4694398 


0.2197783849 
0.5258342731 
0.7940878876 


0.9689372255: 


0.9689372255 
0.7940878876 
0.5258342731 
0.2197783849 


0.2627880673 
0.5824434454 
0.8392844181 
0.9861065906 
0.9452656094 
0.7450396437 
0.4694398388 
0. 1816229357 


0.6384092183 
0.8799441019 
0.9965134344 
0.9154505797 
0.6928832078 


0.4141150246: 


0.1489001176 


0.6928832078 0.7450396437 
0.9154505797 0.9452636094 

1 0.9965134344 0.9861065906 
0.8799441019 0.8392844181 
0. 6384092183 0.5824434454 
0. 3606984983 0.31 0.2627880 
0.122105975 0. 1016466798 


0.08783237415 0.08087246878 
hHPWehEPxw 


Figure 3.8. Calculating the Hamming Window Samples. 


0.8 


σι 


O. 
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o 
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win] 
0.4 
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samples C nl 


Figure 3.9. The Hamming Window. 
and finding the magnitude of the result as shown in a section of the FREQRES 
function (Figure 3.10). The following describes the APL commands which generate 


the frequency responses and Figures 3.11 and 3.12 show the corresponding results. 


HHP<—128 FREQRES hHP 
HHPW—128 FREQRES hHPW 
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VFREQRES(O19 
(01 HWeM FREQRES h;üI0;N;H 


(11 A 
τ] GENERATING FREQUENCY RESPONSE MAGNITUDE AND PHASE FROM TIME DOMAIN 


A 

[31 R COEFFICIENTS 

(41 A 

1 010+0 

(6]  Ne(M-ph)po 
Π.Ν 


9 
C3] HeFFT h 
C9] EWe(M+2)T(XMAGN H) 
(107  PWe(M=2)TC(XPHAS H) 


(113 >O 





(121 A THE FUNCTIONS FFT,XMAGN,XPHASE ARE USED. 

C13] a IN ORDER TO RECEIVE A BETTER RESOLUTION THE COEFFICIENTS ARE ZERO- 
(141 A  PADDET TO 'M' SAMPLES. 

C15] A MeNUMBER OF SAMPLES INCLUDING ZERO PADDING (USE ONLY RADIX TWO NO. 
C16] a heTHE FILTER COEFFICIENTS (VECTOR] 

(171 A HW: FREQUENCY RESPONSE MAGNITUDE (M+2 SAMPLES) 

C18] A PW: FREQUENCY RESPONSE PHASE (M+2 SAMPLES) CRAD. ] 


Figure 3.10. Using the Function FFT. 


2m 2 
di 
m 
2 0.8 
= 
rn 
i 0.6 
T 
d 0.4 
2 
o. 2 
° 





ExP1.064J 


[0] 
Figure 3.11. The Frequency Response Using 51 Coefficients. 


C. EC3410 WORKSPACE 


EC 3410 is a course that gives an introduction to discrete-time random process. 
The topics covered are: description of discrete-time random signals and random 
vectors, linear transformations, sampling of continuous-time random signals, esti- 
mation, and spectral analysis. The EC3410 workspace was developed in order to 


enable the student to solve computer assignments in digital signal processing and 
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VAL +-30p3 





O 20 40 6ο so 
[91 CxPI/“64] 


Figure 3.12. The Frequency Response Using Hamming Window. 


learn to write algorithms in an efficient way without doing excessive programming. 
This workspace contains functions that compute the mean, circular convolution, 
linear convolution, fast Fourier transform, power spectrum and other quantities. 
A full description of the workspace functions is given in Appendix B. 

By using the STATGRAPHICS software package with the EC3410 workspace, 
a graphic capability is achieved. The following is a computer assignment in statis- 
tical dsp and its solution using the EC3410 workspace and the STATGRAPHICS 
software. 

Note that in the actual assignment, students are required to write some of the 
functions that are used below. However developing the functions in APL using the 


tools provided is not difficult. 
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EC 3410 COMPUTER ASSIGNMENT [5] 


. Three data sets SL1.DAT, SL2.DAT, and SL3.DAT are given on a disk. The 
mean of the signal data SL1.DAT is 0.75. Subtract this constant value from 
the signal to generate a new zero-mean signal SLI’ that you will use in this! 
assignment. The other data sets (SL2, SL3) have zero mean. Do the remaining 
parts for each of the data sets SL1', SL2, and SL3. | 


. Using the correlation method, generate a 3 by 3 Toeplitz correlation matrix 
for the signal data. Print this matrix. 


. Solve a set of Normal equations involving the correlation matrix that was just 
generated to obtain the 2nd order linear predictive filter parameters and the 
prediction error variance. | 


. Apply the filter to the original data set and generate the prediction error 
signal. Plot this signal and compute its variance. Does it compare well with 
the theoretical prediction error variance you obtained by solving the normal 
equations? 





Solution 


1. By using the function MEAN, the SL1 data set mean is computed and is 
subtracted from the original SL1. The new data set SL1’ has a zero-mean. 
(Figure 3.13) 


MeMEAN SLI 
M 
2.065552048 


SL1*SL1-M 
MeMEAN SL1 


M 
6.982261991E 17 


Figure 3.13. Using the Function MEAN. 
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ΚΞ SACS SLI 


R1 

8.575461794 8.055195789 7.537100879 
Klecov R1 
K1 


8.575461794 8.055195789 7.537100879 
8.055195789 8.575461794 8.055195789 
7.537100879 8.055195789 8.575461794 


Figure 3.14. Using the Functions SACF to Evaluate the Covariance for 
SL1.DAT. 


R2e3 SACF SL2 


R2 

1.214143223 0.589008547 0.2762704521 
K2€e€COV R2 
K2 


1.214143223 0.589008547 0.2762704521 
0.589008547 1.214143223 0.589008547 
0.2762704521 0.589008547 1.214143 


R3e3 SACF SLI 


R3 

9.410119655 ~8.931635784 8.467433519 
K3ecov R3 
K3 


9.410119655 78.931635784 8.467433519 
~8.931635784 9.410119655 ^78.931635784 
8.467433519 78.931635784 9.410119655 


Figure 3.15. Using the Functions SACF and COV for SL2 and SL3. 


2. By using the function SACF (sample autocorrelation function), the first three 
correlation function lags are calculated and then by using the function COV a 
Toeplitz matrix is generated from the correlation function data. The process 
described in the above paragraph is repeated on the SL2 and SL3 data sets, 
as seen in Figures 3.14 and 3.15. 


3. The function LPFP (linear prediction filter parameters) solves the Normal 
equations and it calculates the 2nd order filter coefficients and the prediction 
error as demonstrated in Figure 3.16. (This could also be done without the 


22 


ale2 LPFP Sil a22 LPFP SL2 


al a2 

1 70.9666951989 0.02913175079 1 “0.4900717314 0.01020142115 
PERR PERR 

1.008111643 0.9283051353 


a3e2 LPFP SL3 


1 0.9593783731 0.01077400633 
0.9325296297 
Figure 3.16. Solving the Normal Equations Using the Function LPFP. 


LPFP function by setting up simple linear equations and solving them as 
described in Chapter II.) 


4. Using the LCV (linear convolution) function, the three data sets are convolved 

with the 2nd order filter parameters that were found in the previous step to 

` provide the prediction error sequence. The prediction errors of the signals 

are determined by computing the initial value of the autocorrelation function 

of the prediction error sequence as seen in Figure 3.17. As can be seen, the 

results are very similar to the prediction error values derived by the analysis 

of Figure 3.16. Figure 3.18 is a plot of the prediction sequence for each data 
set. 


 -— examples of using the EC3410 workspace are provided in Appendix 
D. EC4440 WORKSPACE 

EC 4440 is a course at the Naval Postgraduate School in multidimensional 
digital signal processing. The — deals with the analysis of signals that are 
functions of the two or more independent variables. The course develops both 
time/space and frequency domains approaches and highlights the subjects that 


are different from one-dimensional signal processing. Some of the topics covered 
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ale2 L2FP SLL 


al 
1 70.9666951989 0.02913175079 


PERR 
1.008111643 


a2e2 LPFP SL2 


a2 
L ~0.4900717314 0.01020142115 


PERR 
0.928305135953 


a32 LPFP SL3 


1 0.9593783731 0.01077400633 


PERR 
0293252296297 


Figure 3.17. Calculating The Prediction Error Using The Function PRER. 
are: two-dimensional circular and linear convolution, difference equations, recur- 
sively computable systems, Fourier analysis, and methods of 2-D filter design. The 
EC4440 workspace enables the student to solve computer assignments in multidi- 
mensional signal processing by using functions such as: 2-D convolution, 2-D Fast 
Fourier Transform (FFT) and inverse 2-D FFT, and filter design using the win- 
dowing method or the McClellan transformation. The plots are provided by using, 
as before, the STATGRAPHICS software package. More detailed lists and descrip- 
tions of the workspace functions are provided in Appendices B (DSP library) and 


C (functions listings and samples of runs). 
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PREDICTION ERROR OF SLi 
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Figure 3.18. The Prediction Errors For SL1, SL2, And SL3. 
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The following is an example of how to design a lowpass filter using the workspace 
functions. The design here is done by the McClellan transformation method. 


Again, as in other examples, the students are asked to write some of these functions. 


COMPUTER ASSIGNMENT ON McCLELLAN TRANSFORMATION [6] 





Consider the prototype filter: 


H,(w) 


—T -- ἐπ 4π T ω 


1. Using a 32 point FFT, find the impulse response h(n) and the coefficients 
a(n) in the representation 


N 
it — 3 a(n)T, [cos w] N = 15 


2. Now let F(w,w2) = 5 (1 + cosw; + coswa + cos wa - cosw2) 


Compute H(w,,w2) = ری‎ a(n)T,[F(w,,w2)] and generate both 3-D and 
contour plots of this frequency response. 


Note: You can take advantage of the recursion 
T (Z) = 222 y = =) 


and the built-in recursive nature of APL functions to generate the Chebyshev 
polynomials. 
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Solution 
The function PROFILT is used to specify a 32-sample ideal lowpass prototype 
filter; then by using the function COEFF the impulse response coefficients are 


calculated. As shown in the following steps the COEFF function uses the IFFT 


function to calculate the impulse response, h(n), of the prototype filter. 


VCOEFF(O1V 
[01 a*COEFF Hp;OI0;SIZE:HH:H 
A 


e 
(21 a COMPUTING THE COEFFICIENTS FROM THE PROTOTYPE FILTER (USING IFFT) 


C3] a 
C4] 00410 
(51 SIZE*( eHp)+2 


[61 aeSrzEoo 
| p] 88 HB | 


. The filter coefficients are defined as: 
πω. n=0 
2h(n), n > O 
and are found by performing the following steps in the function COEFF. 
B 


CLO] aC(iSIZEIe€2xHH(LSIZEI 
[111 at01*atol1-2 


(121 >O 
C13] a 
C14] R HpeTHE PROTOTYPE SAMPLES 
[151 a a:THE IMPULSE RESPONSE COEFFICIENTS 
se A THE IFFT FUNCTION HAS BEEN USED. 
17] a 
(181 A Y.KATZIR,I.A.F., OCTOBER 1987 
a 


6191] 


و 


HP©PROTFILT 32 


a¢COEFF HP 


HP 
1 1 1-1 2 1.10700 0 00707010 08070 0 ο ο”. οο τπτ 


a 

0.40625 0.610186859 0.1779849416 “Ὁ. 1365886089 “0. 1508883476 0.0129955742 
0.110335429 0.0464416897 70.0625 ^70.07130573562 0.01466457095 
0.07052675554 0.02588834765 ~0.0504871114 ^70.05298494156 0.01823057754 


Figure 3.19. Calculating the Impulse Response Coefficients. 


PROTOTIPE FILTER H(o) 


UT pe vee Diy LX 
e 
O 


= | 
E 
Qo 





0 1 
Samples Cw] 
Figure 3.20. 1-D Prototype Filter. 


The above calculations are demonstrated in Figure 3.19. (Notice that the prototype 
filter response is defined on the interval 0 to 27 rather than —7 to +r.) Figures 
3.20 and 3.21 are plots of the prototype filter and the prototype filter coefficients. 

The ideal filter would have a step discontinuance at the points w = 0.47 
and 1.67 but this appears in the plot as slanted lines because the plot uses finite 


numbers of samples. 
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Figure 3.21. The Prototype Filter Coefficients. 


The function TRANSFNC generates the transformation function by imple- 


menting the following equation: 
F(w1,w2) = 0.5(-1+ cos wl+cos w2+cos wl-cos w2). 


The surface plot capability of STATGRAPHICS is used to form a 3-D plot anda 
contour plot of the transformation function as seen in Figure 3.22. 

By using the function MCCLEL the frequency response of the desired filter is 
generated. The function generates the nth Chebyshev polynomial, Tn[F(w1, w2)] 
by using the recursion T,(r) = 2(1)T,-1(1) — Tn-2(z). By using the following 


summation: 
N 


>` a(n)T,, [F(w1, w2)| 


nz 
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Figure 3.22. The Transformation Function. 
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the frequency response H(wl,w2) is achieved. The above steps are shown in a 


section of the MCCLEL function. 


VMCCLEL(01V 
(01 Hwwea MCCLEL F;DIO;N;INTER;NUMBER;INTER1;INTERO 
a 


(21 a DESIGNING 2-D FILTER USING MECCLELLAN TRANSFORMATION 
(31 A 

(41 00+0 

(51 INTERO+INTERI©INTER+Hwwe( oF) 00 

C6] NUMBER+*oa 









LOOP:>((N=0),N=1)/A.3 
INTER*e(2xFxXINTER1)-INTERO 

(101 >C 

(111 A: INTERe(oFJol 

(121 INTERI+INTER 

[13] >C 

(141 B:INTEReF 

(151 C:0*'*CALCULATED CHEB. POLYNOM OF ORDER:  , δΝ 
(161 HwweHww*CatN1xINTERJ 
(171 INTERO*FINTER1 
INTER1€INTER 
-(NUMBER»Ne€N*1)/LOOP 












FeTHE TRANSFORMATION MATRIX (CALCULATED USING TRANSFORMATION 


FUNCTION) 
aeTHE PROTOTYPE IMPULSE RESPONSE COEFFICIENTS 


Hww: THE 2-D FREQUENCY RESPONSE. 


Y.KATZIR,I.A.F., OCTOBER 1987 


m 

N 

de 

مه 
DDDDDDDO‏ 


Figure 3.23 shows the frequency response plot of the lowpass filter. In order 


to apply a status report during calculation, a status printout was added to the 


function as seen in Figure 3.24. 


E. DISCUSSION 
The workspaces and functions developed in this thesis provide a complete 
set of tools for solving the problems currently assigned to students in courses EC 


3400, EC 3410, and EC 4400. They also provide capabilities to solve many other 
problems in DSP. 
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Figure 3.23. The Frequency Response of the Lowpass Filter. 
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Figure 3.24. Status Printout During Calculation of Chebyshev Polyno- 
mial. 


Some of the functions developed and included in the workspace are finite 
specialized and were specifically coded as examples of solving the particular DSP 
homework assignments. Others are provided as tools that can be applied to many 
different problems. As already mentioned, it is very easy to expand this set of tools 
to include others. 

In using the workspaces and functions, the instructor will usually provide 
only a subset of the functions developed here to students enrolled in the course. 
Students would then use the functions provided, as well as the general features of 
APL to solve the assigned problems. Since the workspaces are extremely modular, 
an instructor will find it easy to select only those functions that he wants students 


to have and provide others at a later time. 
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IV. SUMMARY AND CONCLUSION 


The purpose of this thesis was to develop a convenient set of software tools for 
students and instructors in order to solve computer assignments in digital signal 
processing on personal computers using APL. 

The main requirement was achieved by developing APL workspaces with sets of 
functions that could be applied to perform DSP operations. Computer assignments 
in EC 3400. EC 3410, and EC 4440 are easily solved by using this software package. 
The user needs to have basic knowledge in APL and needs to know how to use part 
of the STATGRAPHICS package in the APL environment. However, this basic 
knowledge is easily acquired. 

The structure of the APL software is very flexible and convenient to use. By 
combining a few functions, the user can achieve complicated calculations with little 
programming effort. Any workspace can be expanded by the student or professor 
by adding new functions. 

The results of this thesis demonstrate that the APL is easy and convenient 
to use, especially in solving digital signal processing problems. This thesis will 


hopefully make it easier for new students to start using the language. 
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APPENDIX A 
USER MANUAL 


SYSTEM REQUIREMENTS 


The software package contains the following disks: 
. APL*PLUS Executive Program version 5.1 or later by STSC 
. STATGRAPHICS - Statistical Graphics System, version 2.1 or later by STSC 


. Utility workspace, and EC3410, EC4440, and EC3400 workspaces 


STATGRAPHICS requires the following hardware and system software: 

. IBM PC, PC-XT, PC-AT or compatible 

. 512K RAM 

. Keyboard with APL characters 

. Two double-sided, double-density disk drives or one floppy disk and hard disk. 
. Graphics adapter card 

. DOS, version 2.0 or later 

. Dot matrix graphics printer (optional) 


. Math coprocessor (optional but recommended) 


The software package can be operated from floppy disks or from the hard disk. 


B. PREREQUISITES FOR USER 


In order to use the package efficiently, the user should have a basic knowledge 


of APL. 


In addition, it is strongly recommended that a copy of the STATGRAPHICS 


manual [3] be available to the user since this Appendix refers the user to Ref. 3 
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for more explanations. The user should not have loaded any special DOS resident 
programs since this may interfere with the proper execution of STATGRAPHICS 


and cause the system to “hang up.” 


C. HOW TO GET INTO THE DSP LIBRARY 


Follow the instructions below to start the session: 
1. Start the computer using DOS. 


2. If you are using a dual floppy disk drive, place the APL*PLUS executive disk 
in the default disk drive. If you are using a hard disk system, make sure that 
the directory with the APL and STATGRAPHICS program (usually the APL 


directory) is your current directory. 


Enter: APLCOM if you are using APL version 5.1 or 
Enter: A if you are using APL version 6.4. 


These commands call up batch files specifically for each version. A batch file 
contains a set of commands that must be executed before it runs the APL program 
in order to have the software matched with the specific hardware. 


The system will respond as follows: 


VQuQGQ(eS===[[{janha[==*x={K{€V—lC"“***{]_*[_—{]{_{€X[=2T]_TRlHUHHaana={=[=[{=ee —ei<—=_xz=_z_x_zx_z{]_{_§]—>—[]—— Il I> yy oo. a | 
l APL#PLUS PC APPLICATION DEVELOPMENT SYSTEM l 


Version 3.1 Serial Number 134257 Copyright 1995, οἷσο Ine: ۱ 
All rights reserved. Unauthorized reproduction of this software l 
is prohibited and violates U.S. Copyright Laws. APL*PLUS is a i 
registered trademark and service mark of STSC, Inc. i 


ed 
CLEAR WS 11/20/1287 14:31:34 


3. Place the desired DSP library disk in the default disk drive, as for example, 
the disk containing the UTILITY workspace. 


Enter: JLOAD UTILITY 


If you are using a hard disk and want to load the workspace from drive A, 


enter: JLOAD 0 UTILITY. The “0” (zero) stands for drive A, “1” stands for 
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drive B, and “2” stands for drive C (hard disk). The system will respond as 
follows: 


JLOAD 0 UTILITY 
ο UTILITY SAVED EDV LSU 15:12:28 


| NAVAL POSTGRADUATE SCHOOL MONTEREY, CALIFORNIA | 


UTILITIES WORKSPACE 


Y.XATZIR.I.A.F. 
ADVISOR: PROF. C.W. THERRIEN. 
VERSION 1.0 SEPTEMBER 1987 


< s. a an an an an «αρ «πρ «αὐ «πο «παν «αν am am AMD AMD ED ED ED ED ED ED ED A A A A A UD aD UD AUD AUD A A A A UM A A A A aD aD aD a A A m G G «αν «αν «αν e e e 


| | 
| | 
| THESIS Bf l 
۱ l 
| ۱ 
| l 


At this point you can use the UTILITY workspace as described in Appendices 
BSC ۰ 


4. Before terminating the session save your results by entering: 
JSAVE UTILITY or JSAVE 0 WSID (WSID means your file name) 
the system will respond: 


)SAVE 
0 UTILITY SAVED 11/19/1987 00: 24:06 


9. To terminate APL enter: JOFF and you will be automatically returned to 
DOS. 


D. USING STATGRAPHICS WITH APL 


Follow the instructions below to start the session: 


1. Start the computer under DOS. If a graphics printer is attached be sure that 
the printer is set up for graphics printing. (DOS has to run graphics.com when 
system boots.) 


2. Load APL as described in the previous section (C.b.) 


3. After the APL*PLUS PC system is loaded, replace the APL disk with the 
STATGRAPHICS start-up disk and enter: 


JLOAD STATGRAF 
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WARNING: When using the STATGRAF workspace use only the )P COP Y 


command to enter functions or variables into the workspace. 





During the loading procedure, the system will display a copyright banner and 
the following message: 
‘*SYSTEM INITIALIZING. PLEASE BE PATIENT. THIS WILL TAKE A FEW MOMENTS.’ 
When initialization is complete, the system will prompt for questions that are 
related to the computer configuration. Select the appropriate answers (See Ref. 
3, pg. 3-1 to 3-3 for more details.) The system will ask if you want to save the 
settings previously entered in order to have an automatic logon the next time you 


start STATGRAPHICS. The Main Menu will then be displayed as in Figure A.1. 


STATGRAPHICS Statistical Graphics Systemi 





DATA MANAGEMENT AND SYSTEM UTILITIES TIME SERIES PROCEDURES 
A. Data Management L.e Forecasting 
B. System Environment M. Quality Control 
C. Report Writer and Graphics Replay N. Smoothing 
D. Plotter Interface O. Time Series Analysis 
PLOTTING AND DESCRIPTIVE STATISTICS ADVANCED PROCEDURES 
E. Plotting Punctions P. Categorical Data Analysis 
5. Descriptive Methods Q. Multivariate Methods 
G. Estimation and Tasting R. Nonparametric Methods 
H. Distribution Functions S. Sampling 
I. Exploratory Data Analysis T. Experimental Design 
ANOVA AND REGRESSION ANALYSIS MATHEMATICAL AND USER PROCEDURES 
J. Analysis of Variance U. Mathematical Functions 
K. Regression Analysis V. Supplementary Operations 


Use cursor keys to highlight desired section. Then press ENTER. 
lHelp 228 3Savscr 4Prtscr 5 ہ66‎ 7۷353 8Cnd 9Raviaw 10Quit 
INPUT 11/18/87 13:08 STATGRAPHICS Vers. 2.1 ۴ 


Figure A.1. Statgraphics Main Menu. 
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4. When you are in the Main Menu, press the ESC key to enter APL immediate 
execution mode. Execute the function. Enter: RESTART to return to the 
menu-control mode. 


Since the complete set of functions in the STATGRAPHICS system exceeds 
the memory space available on most personal computers, you should load only 
those functions and variables that are necessary to perform the desired tasks. For 


DSP applications the major procedures used are described in the following section. 


D. STATGRAPHICS MAJOR PROCEDURES USED 

After starting STATGRAPHICS, the first display is the Main Menu (See Fig- 
ure A.1.) 

"Ihe Main Menu contains 22 sections under six major categories. 

The Main Menu gives you access to all the data analysis and graphics proce- 
dures. 

A highlighted cursor bar is displayed on the Main Menu screen. In order to 
select any desired Menu Procedure, move the cursor using the arrow keys and press: 
ENTER. STATGRAPHICS procedures are stored in four disks. When using floppy 
disks carefully follow the instruction on the screen. 

The following is a short description of the most frequently used procedures in 
this DSP software package: 

1. System Environments [3] 

The STATGRAPHICS system environment section is under category: 
“DATA MANAGEMENT AND SYSTEM UTILITIES’? 
This section allows you to modify the STATGRAPHICS system environ- 
ment to fit your needs. This section has nine procedures. Those most often used 


are. 
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Screen Options 


= hk کک کک < س‎ ΕΙΕὓΙΕῸῈῸΟΙὗ  ΤΞΓΞςΝΤΤΤΤττοιτιεὗ,ὗὗ ὗὗὗὗ ὗ ὗ ῦ ὗ ὗ ὗ ΣΤ 


Normal Text 


Foreground Color (0-7): 
Background Color (0-7): 
Border Color (0-15): 
Intense (0-1): 

Blinking (0-1): 


Active Text 


Foreground Color (0-7): 
Background Color (0-7): 
Intense (0-1): 
Blinking (0-1): 


OOrrs 


0 
5 
O 
O 


Graphics Mode 


* Resolution (1=Low,2=High): 
Background Color (0-15): 
Palette Number (0-1): 


ron 


* not supported on all boards 


Highlighted Text 


Foreground Color (0-7): 
Background Color (0-7): 
Intense (0-1): 

Blinking (0-1): 


OONO 


Figure A.2. The Screen Options Menu. 


1. Screen option: Allows you to change colors, text, and properties of the graphics 
screen. (See Figure À.2). 


2. Graphic options: Allows you to change sizes, orientation, and positions of the 
graphics screen. (See Figure À.3.) 


Graphics Options 


























Line colors: 123456789ABCDEF12345 
Point colors: 56789ABCDEF1234518838 
Paint colors: 89ABCDEF123451234567 
Text color: 1 
Axis Color Ticmarks Text Size Color Viewpoint 
X 1 In L T Frame: Yes X 600 
4 L In 1 i Grid: Both 1 Y - 1600 
2 l In L L Border: None 1 2 800 
Titlel 1 1 
Title2 l l Graphics Cursor: Yes Graphics: Yes 
Pointsize: .005 Ticmark length: .040 
Window Origin Width 5 ٦ Position: 1 
Horizontal «200 -750 Printing Zoom Factors 
Vertical . 200 . 640 Vertical: 1 Horizontal: 1 


Page Orientation: Vertical 


Figure A.3. Graphics Options Menu. 


2. Plotting Functions [3] 


Plotting functions are found on the Main Menu under the category: 


‘ “PLOTTING AND DESCRIPTIVE STATISTICS’? 
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STATGRAPHICS provides seven plotting procedures you can use to pro- 
duce a variety of plots. The following is an example of using the X-Y Time Plot 
procedure. 


1. Start in the APL environment and define the variables X and Y as follows: 


O STATGRAF APL Ins 
You are now under control of the APL interpreter. 
To resume menu control, type: 
RESTART 
— Xe. 10 
x 
12345257 8 9 10 
— ο 2.5.3 
Ἔ 
Peo 2.5 325 4 5 6 7.5 8.8 9.5 


43 °6 7-5 8.8 9.5 


Enter: RESTART 
Select the X-Y Line and Scatterplots procedures from the Main Menu. 
2. Fill in the Menu entry as demonstrated in Figure A.4. 


(USE CAPITAL LETTERS) 


X-Y Line and Scatterplots 





ΒΕΓ Log 
رن‎ No 
EP y: Y No 


Point labels: 


Point codes: 
Point colors: 
wea. err. X: 
Sea. err. Y: 


Points: Yes 
— Lines: Yes 


Figure A.4. X-Y Plotting Menu 
3. Press the F6 key in order to see the plot. (See Figure A.5.) 
4. Press the F4 key to print the screen. (if graphics printer is attached) 
9. Press the F5 key to make any changes in the plot (including titles). 


6. Press the F10 key to leave the screen and to return to the menu. 
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Figure A.5. The Plot of Y vs X. 


3. Response Surface Plotting [3] 


This procedure produces 3-D surface plots and contour plots, for previ- 
ously generated two dimensional functions. The procedure is under section: “Ex- 


perimental Design” on the main menu under 


“ADVANCED PROCEDURES? ? 
The following is an example of how to use this procedure: 


1. Start with a definition of a matrix of values for a 2-D function in the APL 
environment as in the following example: 


— XX"16 1601256 


A 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 
17 18 19 20 21 22 23 24 25 26 27 ο 111157 
33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 
49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 
65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 
81 82 383 84 85 86 87 38 89 90 91 92 93 94 95 96 
97 98 99 100 101 102 103 104 105 106 107 108 109 110 1i1 112 
113 114 115 116 117 118 119 120 121 1295 τπτ 26 ٣ 
129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 144 
145 146 147 148 149 150 151 152 153 154 155 156 157 158 159 160 
161 162 163 164 165 166 167 168 169 170 171 172 173 174 175 176 
177 178 179 180 181 182 183 184 185 186 187 188 189 190 191 192 
193 194 195 196 197 198 199 200 201 202 203 204 205 206 207 208 
209 210 211 212.213 214 215 216 217 213 219 220 220222 22: 224 
225 226 227 228 229 230 231 232 233 234 235 236 237 238 239 240 
241 242 243 244 245 246 247 248 249 250 251 252 253 254 255 256 
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— 

— 
— 
—— 


2. Restart STATGRAF and select the section “experimental design.” The fol- 
lowing will be displayed. 


| EXPERIMENTAL DESIGN | 


1. Full and Fractional Factorials 
2. Central Composite Designs 
3. Alias Structure 

— > 4. Response Surface Plotting 


3. Select “Response Surface Plotting” and make the entry in the menu as de- 
scribed in Figure A.6. 


Response Surface Plotting 








X-axis minimum: 0 — Y-axis minimum: 0 
X-axis maximum: 16 — Y-axis maximum: 16 
Number of x-axis intervals: 15 ———> Number of y-axis intervals: 15 


Function type: Matrix input 


Figure A.6. Surface Plotting Menu. 


— 


4. Press the F6 key and fill in the entries as shown in the example of Figure A.7. 


Response Surface Plotting 


X-axis minimum: 0 Y-axis minimum: O 
X-axis maximum: 16 Y-axis maximum: 16 
Number of x-axis intervals: 15 Number of y-axis intervals: 15 


Function type: Matrix input 


Grid matrix: XX 


Figure A.7. Selecting the Function Definition. 
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Surface Plotting Options 
— —<—. 2 F“ [@ .nF —əs——s-sy 
του eitle: THE MATRIX XX (EXAMPLE) 
(2 lines) : 
X-Axis Title: X 
Y-Axis Title: Y 
Z-Axis Title: Z 


— > Plot type: Surface 


For surface plots: 


X-axls resolution: 100 X-axis skip factor: 1 
Y-axis resolution: 100 Y-axis skip factor: 1 
Lines parallel to: Both axes Hidden line removal: None 


For contour plot: 
Type: Lines 
Contour levels: 
om 
Number of Divisions: 10 Minimum: 1 Maximum: 256 


Figure A.8. Surface Plot Entry. 


THE MATRIX xx ( EXAMPLE) , 


e 





Figure A.9. Surface 2-D Plot. 


9. Press the F6 key and enter the values for the surface or contour plot as shown 
in the example of Figure A.8. 


6. Press the F6 key to execute the plot, the F4 key to print the plot, and the F5 
key to edit the plot (See Figue A.9). 


T. Press the F10 key to leave this screen and return to the menu. 


WARNING: Since the number of points in 3-D plots may exceed the memory 
space available, do not select hidden or partially hidden lines in complicated 


plots. Select the entry “points” instead of entry “lines” when drawing compli- 
cated contour plots. 
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APPENDIX B 
DIGITAL SIGNAL PROCESSING LIBRARY 
Listed below are functions’ names and brief discussions of the various computer 
algorithms used to accomplish this thesis research. All programs were written in 
APL by the author unless otherwise noted. Program listings and samples of runs 
are given in Appendix C. For each function its name and the page containing its 
listing in Appendix C are given. Then the calling sequence is given, followed by a 


brief description. 


A. UTILITY WORKSPACE 
1. NORD (Pg. 68 
Z+-N NORD S 
This function generates a normally distributed random vector using the 
function UNRD. 


The function inputs are: 
a. N: Number of elements in the vector 
b. S: Vector with two arguments: 
*[S1]: The required mean 


*[S2]: The required standard deviation 


The function’s source is the software package GRAFSTAT. 
2. UNRD (Pg. 68) 


R— N UNRD B 
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This function generates a uniformly distributed random vector. The inputs 
are: 
a. N: Numbers of elements in the vector 


b. B: Required mean 


3. PUTDATA (Pg. 70 
F PUTDATA ’A:\WORK\NAME. EXT? 
This function transfers data (vector, array, etc.) in free format from an 
APL workspace into an IBM-PC DOS file. The function uses APL NATIVE files 
[7]. A bell sounds when the data transfer has been completed. The function inputs 
are: 
a. F: The data to be transferred 


b. ’A:\WORK\NAME.EXT?: The designated native file or path: “A:” Stands for 
the drive type, “WORK” stands for the directory, and “EXT” stands for extension 
of the DOS file. 


WARNING: The user is not permitted to designate a file that already exists 
as a DOS NATIVE file. 


4. GETDATA (Pg. 72) 


FeGETDATA ’A:\WORK\NAME.EXT? 

This function transfers data (vector, array, etc.) from an IBM-PC DOS 
file to an APL workspace. The function uses APL NATIVE files [7]. A bell sounds 
when the data transfer has been completed. 

>A: \WORK\NAME.EXT’: The native file’s name or path has to be transferred: 
“A:” stands for the drive type, “WORK” stands for the directory, and “EXT” stands 


for extension of the DOS file. 
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B. EC3400 WORKSPACE [4] 
1. DIGFREQ (Pg. 75) 


ΤΟ--Ε5 DIGFREQ FC 


This function calculates the digital frequency using the equation: 
η; (B.1) 


f = The analog frequency 
f, = The sampling frequency 
The inputs are: 
a. FS: The sampling frequency [Hz] 
b. FC: The analog frequency [Hz] 
The output 1s: 


TC: The digital frequency [rad.] 


2. LPCOEFF (Pg. 76 
hLPen LPCOEFF TC 
This function calculates the 2I 4-1 'owpass FIR filter coefficients using the 


Frequency Transformation and the equation: 
hr p(n) 2 k/v(n — I)sin([n — I]0:)),n 2 0,1,2,...,2I (B.2) 


where k = 1 and 0, — The digital cutoff frequency 
The calculated coefficients are for a causal nonrecursive filter. In order to 


generate noncausal coefficients use the function SHIFT. The inputs are: 


a. n: Number of required coefficients (Use an odd number.) By choosing a final 
number of coefficients you automatically define a rectangular window. 


b. TC: The digital cutoff frequency [rad.] 
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j. PPCOBPBPPDP (Pe. 77) 
hHPen HPCOEFF TC 


This function calculates the 2J + 1 highpass FIR filter coefficients using 


the Frequency Transformation method and the equation: 
hgp(n) 2 (-1)0*7Phrp(n),n 20,1,2,...,21 (B.3) 


The calculated coefficients are for a causal nonrecursive (linear phase) 
filter. In order to generate noncausal coefficients (zero phase filter) use the function 


SHIFT. The coefficients are normalized. The inputs are: 


a. n: The number of the required coefficients (Use an odd number.) By choosing 
a final number of coefficients you automatically define a rectangular window. 


b. TC: The digital cutoff frequency [rad.] 


4. BPCOEFF (Pg. 78 
hBP~n BPCOEFF TUL 
This function calculates the 27 + 1 bandpass FIR filter coefficients using 


the Frequency Transformation method and the equations: 


Mami e cosa — Γθο]] ΠΜρίπ) n —0,1,2,...,2I (5.4) 
0. = (0, — 04) /2 


0, — (0, + 0.) /2 


0, = The upper digital frequency 

0, = The lower digital frequency 

0. = The digital cutoff frequency 

The calculated coefficients are for a causal nonrecursive filter. In order 
to generate noncausal coefficients use the function SHIFT. The coefficients are 


normalized. The inputs to the function are: 
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a. n: The number of the required coefficients (Use an odd number.) By choosing 
a final number of coefficients you automatically define a rectangular window. 


b. TUL: Vector of two arguments 
TUL [1]: The upper frequency [rad.] 


TUL [2]: The lower frequency [rad.] 


5. BSCOEFF (Pg. 79) 


hBSen BSCOEFF TUL 
This function calculates the 2J + 1 bandstop FIR filter coefficients using 


Frequency Transformation and the equation: 


hgs(0) =1-— hgr(0) (ial 


L2‏ ۲ئ )م وج ط(- > (مہ)ء و ط 
02)/2 — ,0( = .0 
9¢)/2 + ,9( = ,6 
The upper digital frequency‏ = ,0 
The lower digital frequency‏ — ,0 


0. = The digital cutoff frequency 


The calculated nonrecursive coefficients are for a causal nonrecursive filter. 


In order to generate noncausal coefficients use the function SHIFT. ‘The coefficients 


are normalized. The inputs to the function are: 


a. n: The number of the required coefficients (Use an odd number.) By choosing 
a final number of coefficients you automatically define a rectangular window. 


b. TUL: Vector of two arguments 
TUL [1]: The upper frequency [rad.] 


TUL [2]: The lower frequency [rad.] 
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6. Hamming (Pg. 80 
W—HAMMING N 
This function generates an 27 4- 1 point Hamming window using the equa- 


tion: 
"(n — I) 


W(n) = 0.54 + 0.46 cos | T 


| =0,1,2,...521 (B.6) 


The window starts at n = 0 and has 2/ + 1 samples. The input to the function is 
N: The number of samples. 
7. HANNING (Pg. 81) 
W—HANNING N 
This function generates an 2/+1 point Hanning (von Hann) window using 


the equation: 
—I 
W(n) = 0.5 + 0.5 cos [A uns 05:152. 55:24 (B.T) 


The window starts at n = 0 and has 2I + 1 samples. The input to the 
function is: 
N: The number of samples. 
8. BARTLETT (Pg. 82) 
W—BARTLETT N 
This function generates an 21 4- 1 point Bartlett (triangular) window using 


the equation: 
1—f4(n-I) ,0 <n <2I 
W(n) = (B.8) 


0 , otherwise 


The window starts at n = 0 and has 2I + 1 samples. The input to the 
function is: 


N: The number of samples. 
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9. FREORES (Pe. 83) 


HW—¿Z FREQRES h 
This function generates the magnitude and phase of the frequency response 
from the unit impulse response. The functions: FFT, XMAGN, XPHASE are used. In 


order to achieve good resolution the coefficients are zero padded to 'Z' samples. 


WARNING: The ’Z’ needs to be a power of 2. 


The inputs to the function are: 
a. Z: The total number of samples including the zero padding 
b. h: The impulse response coefficients 

The outputs are: 
a. HW: The frequency response magnitude (Z/2 samples) 


b. PW: The frequency response response phase (Z/2 samples) [rad.] 


10. IDEALF (Pg. 84) 


He-TYPE IDEAEF TC 
This function produces the frequency response of various ideal filters. The 
inputs are: 
a. ‘TYPE’: The tvpe of the filter to be generated: 
"LP’: Lowpass 
'HP': Highpass 
'BP': Bandpass 


'BS': Bandstop 


02 


b. TC: A vector containing the following: 
1. In the case of LP or HP: 
TC [1] = The number of samples 
TC [2] = The cutoff digital frequency 
2. In a case of BP or BS: 
TC [1] = The number of samples 
TC [2] = The upper digital frequency 


TC [3] = The lower digital frequency 


The frequencies should be normalized to lie in the interval 0 to 7. [rad.] 
11. FCOEFF (Pg. 86) 
hen FCOEFF HW 
This function calculates the 21 4- 1 point filter coefficients using IFFT (The 
sampling method.) 
The coefficients are for a causal nonrecursive filter. The inputs to the 


function are: 


a. n: The number of the required coefficients (Use an odd number.) By choosing 
a final number of coefficients you automatically define a rectangular window. 


b. HW: The ideal frequency response. 
C. EC3410 WORKSPACE [5, 8] 
1. SACF (Pg. 88) 
R—LAG SACF DAT 
This function computes the sample autocorrelation function (unbiased) 
using the following equation: 


N —1-|k| 
RO =1UN Y^ YugptYa 0Sk<N-1 (B.9) 
n=0 


53 


The inputs are: 
a. DAT: The data block (one-dimensional vector) 


b. LAG: Number of lags to be calculated (the limit is 300 lags) 


2. MEAN (Pg. 90) 


M—MEAN F 
This function calculates the sample mean of a data set using the following 
equation: 


N 
S 7 (B.10) 
= 1 


The input is: 
F: The data block 
3. COV (Pg. 91) 
K—COV R 
This function generates a Toeplitz covariance matrix using covariance func- 
tion. The input 1s: 
R: The covariance function values (vector) 
4. LPFP (Pg. 92) 
at-P LPFP DATA 
This function generates the pth order Linear Predictive Filter parameters 


and calculates the Prediction Error Variance by solving the Normal Equation: 


RO) RO) ... R(p-1) 1 σθ 
R(1) 7s : "I : (B.11) 
ign bams | R(0) ς᾽ 0 


The inputs are: 


a. P: The coefficients order 


94 


b. DATA: The data block 
The outputs are: 

a. a: The filter coefficients 

b. PRER: The estimate prediction error variance 
The following functions are used: 


a. SACF: Calculates the sample autocorrelation of the data block 


b. COV: Generates the Toeplitz matrix of the covariance function 


5. CC (Pg. 93) 


Yex CC h 
This function performs one-dimensional circular convolution between two 


inputs (A(n) (9 z(n)). z and h are the arguments to be convolved. 


WARNING: z and h need to be the same size and one-dimensional, 


otherwise an error message will be generated and the convolution will not 
start. 





6. FCC (Pg. 94) 


Yx FCC h 

This function does one-dimensional fast circular convolution between two 
arguments. The fast operation is established by a shifted matrix multiplication. 
r and h are the arguments to be convolved. Due to memory size this function 
is limited to arguments of less than 80 numbers. If the argument is too long a 


printout message advises the user to use regular circular convolution function. 


WARNING: z and h need to be the same size and one-dimensional, 


otherwise an error message will be generated and the convolution will not 
start. 





55 


Teo τον σι ος] 


Y—zx1 LCV hi 
This function does one dimensional linear convolution between two argu- 
ments by using the following equation: 


OO 


zi(n)*hi(n)= Y zi(m)h(n-m) (B.12) 


m= —oo 
zl and hl are the arguments to be convolved. 
The function uses the CC function (Circular Convolution) and the linear 
convolution is produced by zero padding the arguments. 
8. FLCV (Pg. 96) 
Yexi1 FLCV hi 
This function does one-dimensional linear convolution (same as LCV func- 
tion) by using the function FCC in order to achieve a faster computation time. 
9. FFT (Pg. 97 
Z—FFT X 
This function calculates the one-dimensional Discrete Fourier Transform 
according to the definition 


N —1 
X( = 5 z(n)e GG" eee (B.13) 
zt) 


The input to the function: 


X: Can be either real vector of length N or a two by N matrix (the 1st row 1s the 
real part and the 2nd row is the imaginary part). N must be a power of two. 


The output 1s: 


Z: A complex two by N matrix (real and imaginary) 


The function was written by Professor Paul Penfield, Jr. of M.I.T. 
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10. IFFT (Pg. 98) 


Z—IFFT X 
This function calculates the one-dimensional Inverse Discrete Fourier Trans- 


form according to the definition: 


irem " 
z(n) 2 1/N >` N TED N= (B.14) 


k=0 
The input to the function is: 


X: Can be either a real vector of length N or a two by N matrix (the 1st row is 
the real part and the 2nd row is the imaginary part). N must be a power of 
two. 


The output is: 


Z: A complex two by N matrix (real and imaginary) 


The function was written by Professor Paul Penfield, Jr. of M.I.T. 


11. PSE (Pg. 100) 


INW—ZP PSE DATA 
This function computes the power spectrum estimate (Periodogram) using 


the equation: 
2 


$, (e?) 2 IN(0) — UN x ο (8.15) 





The function uses the FFT APL function and the inputs are: 
a. DATA: The data block (one dimensional vector) 


b. ZP: The length of the zero padding needed for a smooth plot 


WARNING: The DATA plus ZP need to be in length which is a power of two 


Οἵ 


12. PSEB (Pg. 100 
INBW—SEG PSEB DATA 


This function computes the power spectrum estimate with rectangular 


nonoverlapping window using the equation: 


2 





IN") (6) = UN xO (e) 
(B.16) 


M 
B(0) 2 1/M » ٔ ۲ 
=1 
N: The length of the data set 


M: The length of a segment 


Each segment has N' Points: N 2 N' x M. 
This function uses the FFT and the PSE APL functions and the inputs are: 


a. SEG: The length of the segments 
b. DATA: The data block (one-dimensional vector) 


The function automatically zero pads the segment up to 256 points. 


WARNING: In case of error due to memory size, use fewer segments 


13. XSQRT (Pg. 102 


XAZe-XSQRT XR 


This function calculates the complex square root of a complex number. 


XZ = XR1/2 (B.17) 
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14. XTIME (Pg. 102) 


XZ—XL XTIME XR 


This function multiplies two complex numbers 


XZ- XLx XR (B.18) 


15. XMAGN (Pg. 103) 


Z+-XMAGN XR 
This function computes the complex magnitude of a number, the result 1s 


a real number 


Z =|XR| (B.19) 


16. XPHAS (Pg. 103) 


Z—XPHAS XR 
This function computes the phase of the complex number XR, and returns 
a real number. 


Z =Z XR (B.20) 


17. XCONJ (Pg. 104) 


XZ—XCONJ XR 


This function returns the complex congugate of the input. 


18. XEPO (Pg. 104) 


AZ«*—XEPO XR 


This function computes the complex exponential 


XZ = e> (B.21) 
19. PI (Pg. 105) 
Y«—PI 
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This function returns the 7 value. 


20. COS (Pg. 105) 


Y<—C0S X 


This function returns the cosine of the input. 


21. SIN (Pg. 105) 


Y—SIN X 


This function returns the sine of the input. 


D. EC4440 WORKSPACE [9] 
1. CC2D (Pg. 107) 


Y+-X CC2D H 


This function calculates two dimensional circular convolution using the 


equation: 
و‎ = 
Νι--1 ΝΟ --ἱ 
Y (my, nə ) = 3 p? h(ni,n2) z [((n1 — mi))y, ((n2 — ma))y, | (B.22) 


Where subscript Ni and N2 denotes arithmetic modulo Ni and N2. x and h are 


the arguments to be convolved. 





2. FFT2D (Pg. 109) 


ΡΕ ΣΙ 


This function calculates the two-dimensional DFT using the one- dimen- 
sional FFT function. If the input is real, the function converts it to a complex array. 


If the output array is real valued the second plane is eliminated automatically. 


60 


WARNING: The input needs to be a length which is a power of two. 


3. IFFT2D (Pg. 110) 


ITED 
This function calculates the two-dimensional Inverse DFT using one-dimensional 
IFFT function (the same explanation as FFT2D). 
4. SHIFT (Pg. 114) 
Y<SHIFT X 
This function shifts periodic data to an interval starting at the origin. 
5. UNSHIFT (Pg. 114 
Y«—UNSHIFT X 
This function shifts periodic data to an interval centered at the origin. 
6. SHIFT2D (Pg. 115 
Y —SHIFT2D X 
This function shifts two-dimensional periodic data to an interval starting 
at the origin 
7. UNSHIFT2D (Pg. 115) 
Y UNSHIFT2D X 
This function shifts two-dimensional periodic data to an interval centered 
at the origin. Each of the following functions generates a two-dimensional window 


by using a one-dimensional window. Two types of supports are used: 


a. A rectangular region of support which is achieved by computing the outer 
product of two 1-D windows: 


Wr(7n1, 22) = W (ni )W,(n> ) (B.23) 
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b. A circular region of support which is formed by sampling a circularly rotated, 
1-D, continuous function: 


πο. (Vit ent) (B.24) 


8. TEP Paice ea) 
W—R ILPFILT P 


This function generates a circular based rectangular window using the 1-D 


equation: 


W(n) = ٢ eo (B.25) 


, otherwise 


This can also be used for defining the frequency response of an ideal circular sym- 


metrical lowpass filter. The inputs are: 
a. R: The window base radius 


b. P: The dimension of the square matrix in which this window is defined. 


9. HAMMINGR (Pg. 119) 
WRe—II HAMMINGR N 


This function generates a 2-D Hamming rectangular based window using 


the 1-D equation: 


Tae | E" —0.46cos(zn/I) ,|n| « I (B.26) 


, otherwise 


The inputs are: 


a. II: The window base dimension (one-side) 


b. N: The dimension of the square matrix in which this window is defined. 
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10. HAMMINGC (Pg. 121) 


WC«-II HAMMINGC N 
This function generates a 2-D Hamming window with circular base using 


equation (B.26). The inputs are: 
a. II: The window base radius 


b. N: The dimension of the square matrix in which the window is defined. 


11. HANNINGR (Pg. 123) 


WR—II HANNINGR N 
This function generates a 2-D Hanning window with rectangular base using 


the 1-D equation: 


0.5 [1 + οοθ(πη/])] , |n| <I 
0 (B27) 


06 l , otherwise 
The inputs are: 
a. II: The window base dimension (one side) 
b. N: The dimension of the square matrix in which the window is defined. 
12. HANNINGC (Pg. 125) 


WC—II HANNINGC N 


This function generates a 2-D Hanning window with circular base using 
equation (B.27). The inputs are: 


a. II: The window base radius 


b. N: The dimension of the square matrix which this window is defined. 


13. BARTLETTR (Pg. 127 


WR—1I BARTLETTR N 
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This function generates 2-D Bartlett (triangular) window with rectangular 


base using the 1-D equation: 


1 - n Όττο 
win)= distin +0 (B.28) 


0 , otherwise 


The inputs are: 
a. II: The window base radius 


b. N: The dimension of the square matrix in which the window is defined. 


14. BARTLETTC (Pg. 129) 


WC+-II BARTLETTC N 
This function generates a 2-D Bartlett window with rectangular base using 
equation (B.28). The inputs are: 
a. II: The window base radius 


b. N: The dimension of the square matrix in which the window 1s defined. 


15. RECTANR (Pg. 131) 


WR— II RECTANR N 
This function generates 2-D rectangular window with rectangular base 
using equation (B.25). The inputs are: 
a. 11: The window base dimension (one side) 


b. N: The dimension of the square matrix in which the window is defined. 


16. PROTFILT (Pg. 133) 


HPROT—PROTFILT SIZE 
This function produces a frequency response 1-D lowpass prototype filter 
in order to design a 2-D FIR filter using a transformation. The function SHIFT 


has been used. The input 1s: 
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SIZE: The size of the prototype filter (number of samples). 


17. COEFF (Pg. 135) 


A—COEFF HP 


This function finds the coefficients a(n) in the representation: 


N 


H(w) = X` a(n) cos(wn) 


n=0 


where: 


h(o mu 
a(n) = لی‎ n>0 


The function IFFT has been used. The input is: 
HP: The filter’s 1-D prototype. 
18. TRANSFNC (Pg. 137 


Fww+-TRANSFNC SIZE 


(B.29) 


This function produces a transformation matrix to design a FIR filter 


using McCLELLAN transformation. In this example the following function has 


been used: 


F(wy,w2) = 1/2 (—1 + coswı + coswz + cos w cos w2) 


(B.30) 


In order to change F'(w,w2) to a different transformation function, replace 


lines 10, 13, and 14 in the program with the new choices. (See example in Appendix 


C (Pg. 137). The input is: 
SIZE: The number of samples of w, or wa. 


19. CHEB (Pg. 140) 


Y—N CHEB X 
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This function calculates the Chebyshev Polynomials using the following 


recursion: 


T,(2) = 22T,-1(2) = ον (B.31) 


The inputs are: 
a. N: The order of the desired polynomial 


b. X: The argument. 


WARNING: The combination of high order of Chebyshev polynomial and a 


large array to be evaluated will result in a long period of calculating time. 





20. MCCLEL (Pg. 141 
Hww'—a MCCLEL F 
This function computes the 2-D frequency response using McCLELLAN 
transformation according to this equation: 


N 
H(w1,w2) — 3 a(n)T, [F'(Qo1 ,co22 )] (48:92) 


n=0 
In order to indicate the calculating status, the function prompts the present cal- 


culating order of the polynomial. The inputs are: 


a. a: The impulse response coefficients of the prototype filter 





| 
b. F: The transformation matrix (previously calculated using the function TRANSFNC) 
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APPENDIX Ç 
FUNCTION LISTINGS AND EXAMPLES OF USE 





JLOAD 0 UTILITY 
O UTILITY SAVED 12/10/1987 16:37:45 


| NAVAL POSTGRADUATE SCHOOL. MONTEREY, CALIFORNIA | 


| == که سه کے‎ aD €D ab ab ab SO G G G G SP G «πᾶν FB G G G SF G G SF G G SF G G G G s s q “< s. G aD aD G G G G G G G G aD aD A ab ab € ab I 


l UTILITIES WORKSPACE 


l 
l | 
l THESIS BY ۱ 
l Y.KATZIR,I.A.F. | 
l ADVISOR: PROF. C.W. THERRIEN. l 
| VERSION 1.0 SEPTEMBER 1987 | 

JFNS 
GETDATA INFORMATION NORD PUTDATA UNRD 
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co] 
1 
[2] 
[33 
[4] 
[53 
[6] 
E71 
[81] 
[9] 
[10] 
LII) 
[12) 
[13) 
[14] 
[151 
[16] 
[171 
[183 
[151 


[01 
[1] 
C2] 
[3] 
[4) 
ESO 
[6] 
[7] 
C8] 
p9] 
1 


۷80 7۳ 

ZeN NORD P;S;I;T 

A 

A A NORMAL DISTRIBUTION RANDOM VECTOR GENERATOR 
A 


P+2TP,1 

Z©(N) 00 

Iel : 
F10:Te2 UNRD 0.5 0.5 

ΤΕΓΖΧΤ 

5 ۲۷٢3(4 0 

ΣΕΙΟΧΙ551 

Z[I)e€eP[1]J+PL2)JxCT[11x(( s2xgS)+S)*0.5) 
+*F10xt(LIeL1+1)SN 

+0 

A 

A SOURCE: GRAFSTAT 

a NeNUMBER OF ELEMENTS IN THE VECTOR 
8 PeMEAN 

a SeSIGMA (NOT SIGMA*2) 
A 


90 7 
ReN UNRD B 

8 

A A UNIFORM DISTRIBUTION RANDOM VECTOR GENERATOR 
A 

R€CBL11-BL21)*CCN?10000)x2xBE 222410000 

+0 


SOURCE: GRAFSTAT 


A 
A 
A N*NUMBER OF ELEMENTS IN THE VECTOR 
A BeTHE LIMITS (2 NUMBERS) 

a 
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GAUSDe100 NORD O 1 


GAUSD 

1.558274527 ~2.444846636 ~1.09846582 1.122362631 0.580943715 
70.2714227577 0.4143999987 "0.9733729641 “1.021759908 
0.3179468935 1.515647616 0.7495965254 '0.5068186108 0.8851652886 
70.2479112202 ~0.7262416263 ~0.4452313973 "0.6121927891 
70.2086518602 0.5618385258 71.064013245 0.3513946699 1.132830543 
0.1510880125 0.7028398161 ^ 0.05231570544 2.018906718 
0.9240866379 ~1.814329439 0.0351494428 71.806312062 1.028100377 
0.3944470661 0.6391243468 0.8738735133 1.75599265 ~0.3197284973 
70.1366870369 0.6159939661 0.9774365509 ~1.115427909 
70.5496175242 0.03996605725 72.484219992 1.158618634 
71.026625814 1.153535117 ~0.7858857486 0.6349969194 0.8198888762 
7021760331562 0.5627148978  0.127309154 0.5538766995 1.09726933 
70.7311467753 1.404690913 ~0.6203249336 0.2370824053 
71.585608223 ~0.4013458241 ~0.7707101307 ~0.2637646461 
0.9762420248 0.9775308888 1.169573915 0.1594895548 0.5001385714 
1.05500756 0.4509335123 1.271139261 0.8989835792 0.438888127 
71.247429466 0.3241026886 0.3904619568 ~0.4052863418 
0.2929672888 2.565225342 ~0.4576884119 ~1.611439324 ~2.669649901 
70.7594268721 ~0.6750622611 '"1.141739437 2.032501263 
0.9683033719 0.6699827805 0.4199990777 ~2.871065022 1.686296301 
0.02722627758 ~0.9020480504 ~2.053207747 0.08940644014 2.0866291 
0.3651320115 0.8457720835 ~0.1842246182 1.03036774 


UNIFD©100 UNRD 0 100 


UNIFD 

3.22 ~5.78 18.74 "10.04 43.86 91.44 ~63.84 “19.84 “18.58 94.66 ~51.14 
` 46.78 “80.44 723.98 "67.32 69.64 71.06 782.22 “14.14 ”37.6 
766.84 70.92 “93.12 “30.38 18.1 92.06 50.94 “90.96 82.12 24.6 
ΘΟ: Ι.154156 "92.3 61.72 81.84 91.02 “77.52 3.34 “18.2 
26.18 “16.3 “93.02 “75.96 “14.68 66.88 44.58 “75.88 “2.82 “10.68 
787.08 “34.64 81.86 “61.82 “68.26 48.36 35.22 6.46 776.02 “84.16 
26.7 74.7 “62.86 “31.94 0.42 “42.96 96.16 33.74 “52.68 39.3 
"51.16 “8.2 “77.86 "60.26 “3.18 “13.76 16.38 64.04 "5.46 6 
731.52 758.7 “35.88 716.94 “47.26 73.98 “87.34 “94.42 “62.04 
780.84 61.58 “89.98 57.78 “99.9 43.1 60.76 “32.04 “1.38 89.56 
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(0] 

(1] 

(2] 

[3] 

C4] 

ESA 

[6] 

[7] 

[8] 

[9] 

[10] 
ELT 
[12] 
[ΠΠ] 
[14] 
[15] 
[16] 
[177 
C18] 
C19] 
C20] 
[21] 
[22] 
[237 
[24] 
[25] 
C26] 
[27] 
C28] 
C29] 
C30] 
[317 
[32 
[331 
[34] 
L 35] 
C36] 
[E37] 


VPUTDATA(OIJV 
F PUTDATA NAME ;OI0;TIE;FA;NREC;C;FB:FC;BELL:1 
3 


A FUNCTION TO WRITE DATA [F] (VECTOR, ARRAY, ETC. JIN FREE FORMAT 


A INTO AN IBM P.c. D.O.S. FILE ( ' NAME' J 
A 
OLO#Ie1 
ONUNTIE ^71 
TIECAI 
A CREATING A TIED FILE 
NAME ONCREATE TIE 
A WRITING OUT THE SHAPE OF F 
FAe''. 10 0 $oF 
A APPENDING THE SHAPE OF F TO THE FILE 
COTCLF ,OTCNL, FA)ONAPPEND TIE 
A RESHAPING F 
NRECef (x/pF) +4 
Cel*4l(x/CoF2)-1 
Fe (NREC,4)0F 
A WRITING OUT MOST OF THE ARRAY F 
FBes(^1 O JF) 
LOOP: FPBLI;(FBUI;]2'7')/toFBI[I;]]e'-: 
*(NREC>I¢I+1)/LOop 
COTCLF ,OTCNL, FB)ONAPPEND TIE 
A WRITING OUT THE REMAINDER OF F AS THE LAST RECORD 
FCe' ' »SCCTFINREC; 1) 
FCLCFCs'7')/vorFC]e:-| 
(OTCLF,OTCNL, FC)ONAPPEND TEE 
A UNTYING THE FILE 


ONUNTIE 71 
2277-0 7 

BELL 

+0 


۵ 

A THE FUNCTION USES NATIVES APL COMMANDS. 
A 

A Y.KATZIR »L.A.F., 6.1.87. 

A ADVISOR: C.W. THERRIEN. 

A 


` 
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GAUSD PUTDATA 'GAUS.DAT' 


)OFP 


C3\>TYPE GAUS.DAT 


100 
1.558274527 
0.580943715 

“1.021759908 
-0.5068186108 
70.4452313973 
71.064013245 
0.7028398161 
71.314329439 
0.3944470661 
70.3197284973 
71.115427909 
1.158618634 
0.6349969194 
70.127309154 
1.404690913 
70.4013458241 
0. 9775308888 
“1.05500756 
0.438888127 
“0.405288 
“1.611439324 
71.171739437 
0.4199999777 
-9.9020480504 


-2. 4448466136 
“9.271 4227577 
0.3179468935 
0.3851652886 
-0.6121927891 
0.3513946699 
*0.05231570544 
0.0351494428 
0.6391243468 
-0.1366870369 
-0.5496175242 
“1.026564 
0.8198888762 
0.5538766995 
-0.6203249336 
-0.7707101307 
1.169573915 
-0.4509335123 
“1.247429466 
0.2929672888 
“2.669649901 
2.032501263 
“2.8710652 
“2.053207747 


0.3651320115 0.3457720835 -9. 
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“1.09846582 
0.4143999987 
1.515647616 

70.2479112202 

0.2086518602 
1.132830543 
2.018906718 

“1.30631202 
0.87383513 
0.6159939661 


0.03996605725 


1.153535117 
-9.1760331562 
=1. 09726933 

0. 4370824053 
-ο. 4631764646] 

0.1594895548 

1.271139261 

0.3241026886 

2.565225342 
70.7594268721 

0.9683033719 

1.686296301 

0.08940644014 


1.122352631 
70.9733729641 
0.7495965254 
“0.72624163 
0. 5618385258 
0.1510880125 
0.9240866379 
1.028100377 
1. 7553939265 
0.9774365509 
72.484219992 
-9.7858857486 
0.5627148978 
-=0.7311467753 
71.585608223 
0.9762420248 
0.5001385714 
0.8989835792 
0.3904619568 
-0.4576884119 
70.6750622611 
0.6699827805 
0.02722627758 
2.0866291 


1842246182 1.03036774 


VGETDATALO0JV 

FeGETDATA NAME;UIO;DATA;M;SHAP;I;BELL 

A 

A FUNCTION TO READ FREE FORMATED DATA [FJ] (VECTOR,ARRAY,ECT. ) 

A FROM IBM P.C. D.O0.S. FILE ['NAME' ] 

A 
Ὀτοςῖςο 
ONUNTIE "1 

A TYING THE PILE 
NAME ONTIE "1 

A FINDING THE SIZE OF THE FILE. 

MeONSIZE “1 

A READING THE DATA INTO THE WORK SPACE 
FeONREAD 71,82,M 
F[(F=0AV[45])/1p0F1+eD0AV[ 253] 

A EXECUTING THE DATA AS AN APL EXPRESSION 
Fes(“(F=D0AV[10J)v(F=0AV[13]J)v(F=0AV[26J])3)3/F 
ONUNTIE "1 

LOOP: I+I+1 
+(I=10 )/ ERROR 
DATA+IJF 
SHAP+*ITF 
>( (x/SHAP)#(eDATA))/ LOOP 
FeSHAPo DATA 
BELLe-OüTCBEL 
BELL 
30 

ERROR: ‘FILE INCORECTLY FORMATED’ 

-0 


DATA HAVE BEEN CREATED BY APL,PL/I,FORTRAN,EDITOR ETC. 


NUMBERS NEED NOT HAVE DECIMAL POINTS AND ARE DELIMITED BY SPACE. 
FIRST RECORD IS SHAPE OF DATA. 


Y. KALZIR. IAE 6.1.87. 
ADVISOR: C.W. THERRIEN. 


DDDDDDDDD 
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GAUSDeGETDATA 'GAUS.DAT' 


. GAUSD 

1.558274527 72.444846636 ~1.09846582 1.122362631 0.580943715 ~0.27142275 
0.4143999987 ~0.9733729641 ~1.021759908 0.3179468935 1.515647616 
0.7495965254 ~0.5068186108 0.8851652886 ~0.2479112202 ~0.726241626 
70.4452313973 "0.6121927891 ~0.2086518602 0.5618385258 ~1.06401324 
0.3513946699 1.132830543 0.1510880125 0.7028398161 ~0.05231570544 
2.018906718 0.9240866379 “1.814329439 0.0351494428 “1.806312062 
1.028100377 0.3944470661 0.6391243468 0.8738735133 1.75599265 
70.3197284973 ~0.1366870369 0.6159939661 0.9774365509 ~1.115427909 
70.5496175242 0.03996605725 72.484219992 1.158618634 ~1.026625814 
1.153535117 ~0.7858857486 0.6349969194 0.8198888762 ~0.1760331562 
0.5627148978 ~0.127309154 0.5538766995 “1.09726933 "0.7311467753 
1.404690913 ~0.6203249336 0.2370824053 71.585608223 ~0. 4013458241 
70.7707101307 ~0.2637646461 0.9762420248 0.9775308888 1.169573915 
0.1594895548 0.5001385714 1.05500756 "0.4509335123 1.271139261 
0.8989835792 0.438888127 “1.247429466 0.3241026886 0.3904619568 
70.4052863418 0.2929672888 2.565225342 ~0.4576884119 ~1.611439324 
722669649901 ~0.7594268721 ~0.6750622611 “1.171739437 2.032501263 
0.9683033719 0.6699827805 0.4199990777 72.871065022 1.686296301 
0.02722627758 ~0.9020480504 72.053207747 0.08940644014 2.0866291 
0.3651320115 0.8457720835 70.1842246182 1.03036774 


τὸ 


JLOAD O EC3400 
O EC3400 SAVED 12/12/1987 21:23:25 


a A A A OD A OD SD OD کت‎ OD کلت‎ OD تک‎ D SED GP GP Cue aD EDD Q G < هھ تک تک اک تک کہ کت‎ <s ωωα 


| NAVAL POSTGRADUATE SCHOOL MONTEREY , CALIFORNIA | 

۱ هه حص هه هه هه هه هه چ هه هه هه هه هه هه هه هه هه هه هه هه هه هه هه هه هه هه هه هه هه هه هه هه هه هه هه هه هه هه وه هه وه هه هه --ν‏ ئک οὖν s‏ سه وه وه هه هه هه هسه < < | 

| ECE -3400 WORKSPACE ۱ 

۱ ۱ 

۱ THESIS BY l 

۱ Y.KATZIR,I.A.F. ۱ 

۱ ADVISOR: PROF. C.W. THERRIEN. ۱ 

| VERSION 1.0 NOVEMBER 1987 | 

)FNS 

BARTLETT BPCOEFF BSCOEFF CC COS DIGFREQ FCC ΓΕΌΕΙΣ EE 
FLCV FREQRES HAMMING HANNING HPCOEFF IDEALF IFFT INFORMATION Ιον 
LPCOEFF MEAN PI SHIFT SIN UNSHIFT XCONJ XEXPO XMAGN XPHAS 


XSQRT XTIME 


T4 


VDIGFREQ(O1Y 


(01 TC*FS DIGFREQ FC 
C13] a 
(21 a CALCULATING THE DIGITAL FREQUENCY 
(31 A 
C4] TCeC2xPIxFCO*FS 
(51 20 
(61 A 
(71 a FSeSAMPLING FREQUENCY (H2] 
(8] A FCeTHE CORRESPONDING ANALOG FREQUNCY [Hz] 
(91 A TC:THE DIGITAL FREQUENCY (RAD. 1 
(101 ñ 
(111 A Y.KATZIR,I.A.F. NOVEMBER 1987 
(121 A 
FC*10000 
FSe50000 
TCeFS DIGFREQ FC 
TC 
1.256637061 
TC+PI 
0.4 


Τὸ 


(0j 
(11 
(21 
(31 
C4] 
(51 
(61 
(71 
(81 
(9] 
(10] 
(11] 
(121 
[13] 
C14] 
C15] 
C16] 
[171 
(181 
t191 
C20] 
C21] 


VLPCOEFF(OIV 


hLPen LPCOEFF TC:I;0I0 


DDD BD 


DDDDDDDD 


CALCULATING THE ‘n’ LOWPASS FIR FILTERS COEFFICIENTS USING FOURIER 
METHOD 


020+7 

nen-1 

nel (n+2) 

hLPenoed 

hLP*e(1+(PIx(1n)))xSIN((1n)3xTC) 

ALP*(TC+PI),hLP 

hLP+*(nt(onLP))3),hL2 

üroeo 

+0 
THE COEFFICIENTS ARE FOR CAUSAL NONRECURSIVE FILTER. 
THE COEFFICIENTS ARE NORMALIZED. 
neoDD NUMBER OP REQUIRED COEFFICIENTS (RECTANGULAR WINDOW) 
TCeDIGITAL CUTOFF FREQUENCY (CRAD.] 


Y.KATZIR.I.A.F. NOVEMBER 1987 


۲6-047 
ne21 


hLPen LPCOEFP TC 


hLP 
71.559907461£^"17 70.0336367435 ~0.02338723209 0.02672826525 0.05045511524 


71.559907461E"17 70.07568267286 70.06236595225 0.09354892838 
0.3027306915 0.4 0.3027306915 0.09354892838 ~0.06236595225 
70.07568267286 71.559907461E 17 0.05045511524 0.02672826525 
70.02338723209 70.0336367435 71.559907461E 17 


TO 


[01 
(11 
C2] 
C3] 
(41 
(51 
(61 
(71 
(8] 
[91 
(101 
ELI 
[121 
[13] 
(141 
[151 
C16] 
C17] 
[181 


YVHPCOEFF(U1V 
hHPen HPCOEFF TC;UIO;N;hLP 


CALCULATING THE ’n’ HIGHPASS FIR FILTER COEFFICIENTS USING FOURIER 
METHOD 


ο ο ο Ὁ 


ατοςο 

Nen 

TCePI-(TC) 

hLPen LPCOEFF TC 
hHPe( C-1)*tN)xhLP 
+0 


THE COEFFICIENTS ARE FOR CAUSAL NONRECURSIVE FILTER. 

THE MAGNITUDE IS NORMALIZED. 

neODD NUMBER OF REQUIRED COEFFICIENTS (RECTANGULAR WINDOW) 
TC+DIGITAL CUTOFF FREQUENCY CRAD. ] 


Y.KATZIR,I.A.P.,NOVEMBER 1987 


DDDDDDDOD 


TCe0.4xPI 
ne21 


hHPen HPCOEFF TC 


hHP 


72.339861191E"17 0.0336367435 0.02338723209 70.02672826525 ~0.05045511524 


72.339861191E" 17 0.07568267286 0.06236595225 ~0.09354892838 
70.3027306915 0.6 70.3027306915 70.09354892838 0.06236595225 
0.07568267286 72.339861191E"17 70.05045511524 70.02672826525 
0.02338723209 0.0336367435 72.339861191E 17 


TT 


to] 
(11 
(21 
(31 
(41 
[5] 
C6] 
(71 
(8] 
(9] 
(10] 
t111 
(121 
(13] 
(14] 
(151 
C16] 
(171 
(131 
C19] 
C20) 
C21] 
(221 


VBPCOEFF(C1V 


hBPen BPCOEFF TUL;UIO;TC;TO;hLP 


Ὁὸ Ὁ Ὁ Ὁ ο ο Ὁ‏ مھ 


CALCULATING THE ‘n’ BANDPASS FIR FILTER COEFFICIENTS USING FOURIER 
METHOD 


00+10 

TCeCTULLO1-TULCL11)*2 
TO+(TUL([01+TUL[1])+2 

hLPen LPCOEFF TC 

nen-1 

net (n+2) 

hLPenihLP 
ABP*(2xC0S((1n+1)xT0))xhLP 
hBPeCntCOnBP)) ,.RBP 

+0 


THE COEFFICIENTS ARE FOR CAUSAL NONRECURSIVE FILTER. 

THE COEFFICIENTS ARE NORMALIZED. 

neOoDD NUMBER OF REQUIRED COEFFICIENTS (RECTANGULAR WINDOW) 
TULeDIGITAL UPPER AND LOWER CUTOFF FREQUENCIES (TWO ELEMENTS) [RAD] 


Y.KATZIR,I.A.P. NOVEMBER 1987 


TULe200 
TUL(11€*0.6xPI 
TUL(21€0.4xPI 
ne21 


hBPen BPCOEFF TUL 


hBP 


77.799537305E"18 1.204383817E "17 0.04677446419 73.152585957E"17 70.1009102 


3.895627305E^17 0.1513653457 73.156840459£717 ~0.1870978568 
1.207227699E"17 0.2 1.207227699E ^17 ^"0.1870978568 73.156840459E 17 
0.1513653457 3.895627305E717 ^70.1009102305 73.152585957E 17 
0.04677446419 1.204383817E"17 77.799537305E"18 


78 


(01 
(11 
(21 
(31 
(41 
(51 
C6] 
(71 
(81 
(91 
C10] 
(111 
(121 
E13J 
(141 
C15] 
C16] 
C17] 


۹8 0 ۳7 


hBSen BSCOEFF TUL;OIO;hBP 


CALCULATING THE 'n' BANDSTOP FIR FILTER COEFFICIENTS USING FOURIER 
METHOD 


πτοςο ° 
hBPen BPCOEFF TUL 
hBSe-hnBP 
nen-1 
hBS(n*21€1-hBPCn*2] 
30 
THE COEFFICIENTS ARE FOR CAUSAL NONRECURSIVE FILTER. 
THE COEFFICIENTS ARE NORMALIZED. 
ر00 ٭ہ‎ NUMBER OF REQUIRED COEFFICIENTS (RECTANGULAR WINDOW) 
TULeDIGITAL UPPER AND LOWER CUTOFF FREQUENCIES (TWO ELEMENTS) (RAD) 


Y.KATZIR,I.A.PF.,NOVEMBER 1987 


TUL(11€0.4xPI 
TUL(21€0.2xPI 
ne21 


hBSen BSCOEFF TUL 


hBS 


7.799537305E"18 0.01284809274 70.01445410434 ”0.06997550689 70.08163809137 


2.341241641E"17 0.1224571371 0.1632761827 0.05781641735 ~0.115632834 
0.8 70.1156328347 0.05781641735 0.1632761827 0.1224571371 
2.341241641E"17 “Ὁ. 08163809137 70.06997550689 70.01445410434 
0.01284809274 7.799537305E”18 


19 


(0] 
(1] 
C2] 
C3] 
C4] 
C5] 
(61 
(71 
(31 
(91 
(101 
(113 
C12] 
[13] 
C14] 
C15] 
C16] 
(171 
C13] 
C19] 
C20] 
C21] 


VHAMMING( O17 
WeHAMMING N;M;üIO 

A 

A GENERATING THE 'N' WIDE HAMMING WINDOW 
A 

στου 

WeNoO 

3((N+2)A4(Ff (N+2)))/0DD 

A EVEN NUMBER OF SAMPLES IN WINDOW 
Me (N+2) 

We0.54+0.,46xC0S( ((PIx2)x1M)+N) 
We(oWw) , WwW 

0 
ODD: Me(N+1)4+2 
We0.54+0.46xC0S(((PIx2)x1M)+N) 
We(eW) ,1LIW 

0 


THE WINDOW STARTS AT t=0 AND HAS 'N' SAMPLES. 
NeNUMBER OF SAMPLES [VECTOR] 


Y.KATZIR,I.A.F..NOVEMBER 1987 


DDDDDO 


Ne31 


WeHAMMING N 


W 


0.083236011124 0.101095942721 0.1378005566 0.1909712636 0.258431248 0.3374186 


0.4246998352 0.5167013823 0.6096567777 0.6997604163 0.7833234448 
0.8569247828 0.917551183 0.9627205933 0.990583773 1 0.990583773 
۵۰99527205933 ۵۰917551183 ١٠٥856924782383 0.7833234448 0.6997604163 
0.6096567777 0.5167013823 0.4246998352 0.3374186903 0.258431248 
0.1909712636 0.1378005566 0.1010959421 0.08236011124 


80 


VHANNINGCOIV 
(01 WeHANNING N:M;QIO 


(11 A 

[2] a GENERATING THE 'N' WIDE HANNING WINDOW 
(31 A 

(4] 010+0 

[5] WeNod 


(61 9CCN*20 Cl CN*222)/0DD 

(71 A EVEN NUMBER OF SAMPLES IN WINDOW 
[81 Me(N+2) 

(9] We0.5+0.5xCO0S(((PIx2)x1M)+N) 

(101 We(eWw),W 

(111 >O 

C12] ODD: Me(N+1)+2 

C13] We0.5+0.5xC0S(((PIx2)x1M)+N) 

C14] WeCcoW),11W 

(151 >O 

(161 
C171 
C18] 
C19] 
(201 
C21] 


THE WINDOW STARTS AT t=0 AND HAS 'N' SAMPLES. 
NeNUMBER OF SAMPLES (VECTOR] 


Y.KATZIR,I.A.P.,NOVEMBER 1987 


DDDDIDOD 


Ne31 
WeHANNING N 


wq 

2.565338304E 3 0.0229303718 0.06282669193 0.1206209387 0.1939470087 
0.2798029242 0.3746737339 0.4746754156 0.5757138888 0.6736526264 
0.7644820052 0.8444834595 0.9103817206 0.9594789058 0.9897649706 1 
0.9897649706 0.9594789058 0.9103817206 0.8444834595 0.7644820052 
0.6736526264 0.5757138888 0.4746754156 0.3746737339 0.2798029242 
0.193939470087 0.1206209387 0.06282669193 0.0229303718 2.565338304E 3 


81 


VBARTLETT(ü1V 
(01 We€BARTLETT N;M;UIO 


(1] a 

[2] A GENERATING THE ’Ν’ WIDE BARTLETT (TRIANGULAR) WINDOW 
[3] Rñ 

C4] 004-0 

[5] WeNoO 


(61 92CCN*2)3CT CN*2) ) )/ODD 

(7] A EVEN NUMBER OF SAMPLES IN WINDOW 
(31 Me(N+2) 

(9] We(eo(1-C(1M)+M)),1-(1M)+M 

[101 40 

(111 ODD: Me(N+1)+2 

[121 Wel-(1M)+M 

C13] WeCeW) .1LW 

(141 >O 

C15] 
C16] 
[17] 
C13] 
C19] 
(201 


THE WINDOW STARTS AT t:0 AND HAS 'N' SAMPLES. 
N©NUMBER OF SAMPLES [VECTOR] 


¥.KATZIR.I.A.F. NOVEMBER 1987 


>2 p pp pp p 2 


Ne31 
WeBARTLETT N 


W 

0.0625 0.125 0.1875 0.25 0.3125 0.375 0.4375 0.5 0.5625 0.625 0.6875 0.75 
0.8125 0.875 0.9375 1 0.9375 0.875 0.3125 0.75 0.6375 0.625 0.5625 Cf 
0.4375 0.375 0.3125 0.25 0.1875 0.125 0.0625 


82 


VFREQRES COV 
(01 HWeM FREQRES h;ü1I0O;N;H 


(11 A 

[2] A GENERATING FREQUENCY RESPONSE MAGNITUDE AND PHASE FROM TIME DOMAIN 
C3] a COEFFICIENTS 

[41 a 

(51 nroeo 

(61 N*(M-ph)p0 

(7] heh, N 


C8] HeFFT h 
C9] HWe(M+2)T(XMAGN H) 
[101 PWe(M+2)?(XPHAS H) 


[111 70 

C12] A THE FUNCTIONS FFT,XMAGN,XPHASE ARE USED. 

C13] a IN ORDER TO RECEIVE A BETTER RESOLUTION THE COEFFICIENTS ARE ZERO- 
[14] A PADDET TO 'M' SAMPLES. 

[15] R MeNUMBER OP SAMPLES INCLUDING ZERO PADDING (USE ONLY RADIX TWO NO. 
C16] A heTHE FILTER COEFFICIENTS [VECTOR] 

C17] a HW: FREQUENCY RESPONSE MAGNITUDE (M*2 SAMPLES) 

C18] a PW: FREQUENCY RESPONSE PHASE (M*2 SAMPLES) (RAD.]1] 

C19] a 

(201 A Y.KATZIR,I.A.F,NOVEMBER 1987 

(21] A 


hLPe21 LPCOEFF (0.4xPI) 


HLP-128 FREQRES hLP 


HLP 

0.9567807992 0.9619585014 0.9763182673 0.9965906571 1.018111173 1.03582952 
1.045421125 1.044271236 1.032124069 1.01125183 0.986091438 0.9624009 
0.9460845406 0.9419036347 0.9523233652 0.9767236974 1.011145443 
1.048645824 1.080227975 1.096202394 1.087755987 1.048461283 0.975462 
0.8701304604 0.7380554233 0.5883839646 0.4325882327 0.282876649 
0.1505027492 0.04424706035 0.03068671843 0.0731944037 0.08614789716 
0.07567134118 0.04996228114 0.01787457278 0.01248921845 0.0349417940 
0.04595168274 0.04486445606 0.03360579506 0.01595981635 3.413311717E 3 
'0.02001949514 0.03043376941 0.03290460379 0.0275668884 0.01625216492 
1.963679943E 3 0.01185437487 0.02209076454 0.02662206686 0.024719091 
0.01713279563 5.857126075E73 6.366569896E^3 0.01668743018 0.02278343 
0.02336291733 0.01842504363 9.228990525E^3 2.014864013E 3 0.01265918 
0.02022375085 


83 


(91 
(17 
C2] 
[31 
C4] 
CS] 
(61 
(71 
(8] 
(9] 
(10] 
(111 
C12] 
[13] 
C14] 
C15] 
C16] 
[177 
(18] 
C19] 
C20] 
C21] 
C22] 
(231 
C24] 
C25] 
C26] 
C27) 
C23] 
C29] 
C30] 


VIDEALF(OIJV 
HeTYPE IDEALF TC;DIO 
A 
A GENERATING IDEAL FREQUENCY RESPONSE 
A 

Qroeo 

87+6 0 

*C((x/CTYPE=’LP’ ))v(x/C(TYPE=’HP’ )))/LP 
-COx/CTYPES ' BP' )) v Cx/ CTYPE* ' BS' )))/BP 

'IMPROPER ARGUMENT,USE ONLY: LP, HP , BP ,OR BS. 
+0 
LP: HCuCh CCTCOOIJ+2)xTCO1I)+1LIJeL 

He(TCCOJ+2)TH 

»(x/C(TYPE="LP* ))/CONTINUE 

He~H 

CONTINUE 
BP:HC CUCL CC TCCOIJF+2)xXCTCOLI-TCCZIIIIIMCL CC CTCCOIJ+2)xTCLC2II+1) I IJeL 
HeCTCCOJ+2) 7TH 

*(x/(TYPE=’BP’ ))/CONTINUE 


He-H 

CONTINUE: HeH, OH 

+0 

A 

a ‘TYPE’ ©THE TYPE OF FILTER TO BE GENERATED (LP.HP,BP,8S] 
A TCCOJ*NUMBER OF SAMPLES (USE ONLY RADIX TWO NO.) 

A TC(1]+THE CUTOFF DIGITAL FREQUENCY FOR LP OR HP PILTER OR 
A THE UPPER DIGITAL FREQUENCY FOR BP OR BS FILTERS 

A TCC2]*THE LOWER DIGITAL FREQUENCY FOR BP OR BS FILTERS 
A THE FREQUENCIES SHOULD BE A FRACTION OF PI CRAD.]. 

A 

A Y.KATZIR ,I.A.P. NOVEMBER 1987 

A 


34 


IDEAL? TC 


76) 158 
TC(21€0.6 


He’ HP’ 


οσο 
οσο 
O + O 
O HO 
Ordo 
O 10 
O HO 
O HO 
O + O 
Oreo 
Ο 
Ono 
O + O 
O HO 
Odo 
O HO 
O HO 
O 10 
Ono 
ہ ہہ‎ 
ںی ہپ "ہو‎ 
پت و بہ "9و‎ 
O HMH O 
O HHO 
O HA O 
O HH O 
O HFA O 
O HMH O 
Ο "0 
O HMH O 
O MH O 
O HAO 
O + 10 


MH O HH O 


IDEALF TC 


TC#300 
τος 119128 
τοί 2190. 3 
TCL31*0.12 
He’ BP’ 


ο ο A 
ο ο η 
O O O 
O OO 
O O O 
O O O 
O O O 
O O O 
OOO 
OOO 
OOO 
OOO 
OOO 
O O O 
O O O 
O OO 
O O O 
O O O 
coo 
eoo 
1000 
mH OOO 
1000 
ہپ‎ OOO 
1000 
1000 
1000 
^ O O O 
HOO 
O O O + 
O OO 
O O ου 
O O O + 


i O O O A 


85 


(01 
[11 
C2] 
[5η 
C4] 
[5] 
C6] 
(71 
(81 
(91 
(10] 
(11] 
(127 
[13] 
(14] 
(15] 
C16] 
(171 
(181 
C19] 


VFCOEFF(OI]V 
hen FCOEFF HW;Ü0IO;H;m;SIZE 
A 

a CALCULATING THE 'n' FILTER COEFFICIENTS USING IFFT 
a 

00-0 

nen-1 

SIZEet (n+2) 

meSIZE*1 

HeIFFT HW 

heC1,m)TH 

hemoh 

heCSIZETCOh2)0,h 

«0 


THE COEFFICIENTS ARE FOR CAUSAL NONRECURSIVE FILTER. 
neODD NUMBER OF REQUIRED COEFFICIENTS (RECTANGULAR WINDOW) 
HWe THE IDEAL FILTER SAMPLES (RADIX TWO) 


Y.KATZIR,.I.A.F.,NOVEMBER 1987 


H 
0000011111111111000000000000000000 
00000000000000000000000 0000000000 O 
0000000000000000000000000000000011 
1111 117070 *0 "05050 0/9 

ne21 


hen FCOEFF H 


h 


0.02511120991 0.04315473177 0.03628637547 76.290075645E^3 ~0.07065358471 


70.1217862264 70.1233903877 70.06216894475 0.03973030337 0.133830098 
0.171875 0.1338300989 0.03973030337 70.06216894475 70.1233903877 
70.1217862264 70.07065358471 76.290075645E^3 0.03628637547 0.0431547 
0.02511120991 
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JLOAD O EC3410 
O EC3410 SAVED 12/11/1987 12:48:01 


ECE -3410 WORESPACE 


Y. KATZIR,I.A.PF. 
ADVISOR: PROF. C.W. THERRIEN. 


۱ 
۱ 
I THESIS BY 
۱ 
I 
۱ 


VERSION 1.0 SEPTEMBER 1987 
JFNS 
cc Cos COV FCC FET FLCV IFFT INFORMATION LCY 
LPPP MEAN PI PSE PSEB SACF SIN XCONJ XEXPO XMAGN 


XPHAS XSQRT ALIME 


87 


(01 


(16) 


VSACF(O1V 
ReLAG SACF DAT ;:010;MEAN;N 
H 
A SAMPLE AUTOCORRELATION FUNCTION (UNBIASED) 
A 
ατοςξςο 
ΝΕρΌΑΤ 
ReLAGo0 
DAT+DAT,R 
LOOP: R(K]+(+/DATx(-K) ODATIJ=N 
+( LAG>K+K+1) / LOOP 
30 


Y.KATZIR,I.A.F. JULY 1987 


DAT+DATA VECTOR (ONE DIMENSIONAL) 
LAGNO. OF LAGS (MAX. 300) 


DDDDDD 


Rxx+100 SACF GAUSD 


Rxx 


1.1816 70.011604 “0.066841 “0.14886 “0.11086 70.050135 70.28244 


0.12983 0.011973 72.0869E"3 "0.052813 "70.066688 0.1414 0.16079 
0.021921 "0.050052 "0.03314 70.12319 "0.11111 "0.087361 
7.9253E°3 0.047914 0.082207 0.022012 70.072971 0.047907 
70.028347 0.03337 "70.091018 0.11214 4.3453E 3 "0.020604 0.08519 
70.13963 0.038409 "0.059126 70.023359 "0.04759 0.03288 0.23321 
70.021454 0.051991 70.049765 "0.025113 0.059572 0.11704 0.11551 
70.11528 "0.058985 "0.038291 "0.077532 0.11274 0.10221 0.18426 
70.09541 70.14171 8.3982573 "0.10611 "0.13464 0.037944 0.035897 
0.13416 70.081343 "0.028463 0.13128 0.026911 0.011188 "0.048038 
0.048817 "0.023835 "0.14964 "70.071093 0.084411 0.04223 0.055516 
0.08221 70.011808 "0.03549 '"4.7687E 3 "0.011637 0.010587 
7.3499Ε”-3 0.033257 70.065615 70.036728 0.038447 0.014501 
0.029666 0.042598 "0.096824 5.4581E"3 0.052421 0.047676 


70.044033 "0.049402 0.016731 71.8189E'"3 6.5336ΣΕ 3 "0.028012 
0.015022 
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CORRELATION OF HITZ HOISZ. 
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(01 
C17 
C27 
[3] 
C4] 
(51 
[5] 
[7] 


VMEANCOJV 
MeMEAN F . 

A CALCULATE THE MEAN 
Me(+/+/F)+x/pF 


A 
A SOURCE:MDSP W. SPACE 
A FeDATA VECTOR 
A M:THE MEAN (NUMBER) 
A 

24+10 

X 

45678910 11712 13 
28 29 30 31 32 33 34 35 
50 51 52 53 54 55 56 57 
72973 74 ٦ 77 78 79 
94 95 96 97 98 99 100 
XXeMEAN X 

XX 

Xex-xx 

X 

748.5 747.5 46.5 745.5 
2238.5 175357 5 535 6 735.5 
728.5 275 726.5 725.5 
1855 αρ στις ς “15.5 
7.5 76.5 5 5 “4.5 73.5 


32.5 33.5 34.5 35.5 36.5 37.5 38.5 39.5 40 


XXeMEAN X 
XX 


l4 15 
36 37 
58 59 
80 81 


90 


16 17 
38 39 
60 61 
82 83 


18 19 
40 41 
62 63 
84 85 


24 25 26 27 
46 47 48 49 
68 69 70 71 
90 91 92 93 


22 23 
44 45 
66 67 
88 89 


20 21 
42 43 
64 65 
86 87 


-39.5 
29.5 


740.5 
-30.5 
"20.5 719.5 

"10.5 “9.5 B5 
2-5 3.5 4.5 5.5 


14.5 15.5 16.5 17.5 1875 


-5 28.5 29.5 30.5 31.5 
¿5 41.5 42.5 43.5 44.5 


vCOVLDJV 
[OMM KeCOV R;UION:;I 


[1] A 

(21 A GENERATE TOEPLITZ COVARIANCE MATRIX FROM COVARIANCE FUNCTION 
(31 A 

[41 üIoeieo 

[5] NeoR 


C6) Ke (N,N) 00 

[71 Re€CO11R),R - 

[8J LOOP: K(IsJeCeN)TC-IIIR 
[9] »>(N>I*I+1)/LOOP 


[10] A 

C113 A SOURCE: MDSP W.SPACE 

[121] A ReTHE CORRELATION FUNCTION (VECTOR) 
C13] a 


RCL1€3 SACF CL1 


RELL 
7.3328 6.8337 6.4035 


CMCL1<4COV RCL1 


GMCL1 
7.3328 6.8337 6.4035 
6 8337 7.3328 6.8337 
6.4035 6.8337 7.3328 
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CO) 
DI 
(2) 
[37 
(41 
(539 
(61 
C71 
[8] 
[9] 
[10] 
LIII 
C12] 
(137 
(14] 
(15] 
(16] 
[17] 
(18] 
C19] 
C20) 
[21] 
[221] 
[23] 
[24] 


VLPFP(O]V 
atP LPFP DATA;UIO;R;K;MDATA;M 

a 

A GENERATE THE Pth ORDER LINEAR PREDICTION FILTER PARAMETERS 
A AND THE PREDICTION ERROR VARIANCE 

A 

OL0+0 

a+(P+1)p0 

MeMEAN DATA 

MDATA€DATA-M 

Re(P+1)SACF MDATA 

KeCOV R 

Δ61,Δ6-(14Β)Β(Ρ Ρ)ὺΚ 

PERRe+/(axR) 

30 


THE FUNCTION SOLVES NORMAL EQUATIONS USING THE SACF,COV,AND MEAN 
FUNCTIONS. 


Y.KATZIR,I.A.F.,AUGUST 87. 


a: FILTER COEFFICIENTS 

PERR: THE ESTIMATE PREDICTION ERROR VARIANCE 
PeTHE FILTER ORDER 

DATA+THE DATA VECTOR 


DDDDDDDDDDOD 


8-3 LPFP CL1 


a 
1 "0.897346235 -70.0173086659 ~0.02113346634 


PERR 


0.9625826774 
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vcccau1V 
ھ 2:7017 22۰۰770060 پر ے۴۱‎ 260 

3 

A 1-D CIRCULAR CONVOLUTION 

2 

020-10 

SNN+p(NN+*pH) 

NepX 

> (NANN)/ERROR 

+(SNN #1) / ERROR 

YeNoO 
LOOP:Y[IJ1€e€*/Xx(-1*I0900H 

-(N2Ie€1*1)/LOOP 

+0 
ERROR: 'ARGUMENTS ARE NOT THE SAME SIZE OR 1-D. 
+0 

a 

A Y.KATZIR,I.A.F,AUGUST 1987 

8 

A ARGUMENTS NEED TO BE THE SAME SIZE AND 1-D. 
A X AND HeARGUMENTS TO BE CONVOLVED (VECTOR) 
ñ 


YeX CC H 


E 
10 


LONG 
'4 5 6 7 8-9 10 


X'CC LONG 


ARGUMENTS ARE NOT THE SAME SIZE OR 1-D. 
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(01 


۹", 7 
YeX FCC H;UIO;N;NN;SNN;A 


[118 
[21 aA FAST 1-D CIRCULAR CONVOLUTION 
(31 A 
C4] 00+0 
[5] SNNe€o CNN*oH) 
[61 NeoX 
C7] > (NN>80)/STATM 
(8] »(N#NN) /ERROR 
[97] -CSNN?1)/ERROR 
[101 A<(N,N)o*$H 
[111 ΑΕ(-ιΝ}ΦΆ 
[12] Y€19CA*.xX) 
[131 +0 
[141 ERROR: 'ARGUMENT ARE NOT THE SAME SIZE OR 1-D.’ 
[151 +0 
C16] STATM: ’PLEASE USE CC OR LCV FUNCTIONS.’ 
(171 >O 
C18] a 
C19] A Y.KATZIR.I.A.F.,AUGUST 1987 
(201 A 
(211 A ARGUMENTS NEED TO BE THE SAME SIZE AND 1-D. 
(221 A LENGTH OF ARGUMENT VECTOR NEED TO BE SHORT. 
(231 A X AND HeARGUMENT TO BE CONVOLVED 
C24] a 
x 
3 2 
H 
3 2 
Zex FCC H 
2 = 
23939 390 
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a —r.. 


VLCVL(DO0JV 


ο. Yexl LCV Hl1;U0IO;N1;N2:;N;X2;H2;X;H 
(11 Rñ 
(21 A 1-D LINEAR CONVOLUTION 
ESO A 
C4] 00-0 
c5] NlepX1 
(61 Ν2ερΗι 
[73] NeN1+N2-1 
(81 XZeNTX1 
(91 82-1 
(101 YeXZ CC HZ 
Cll] «0 
[121 ñ 
(131 A LINEAR CONVOLUTION USING CC FUNCTION 
[14] A 
C15] a Y.KATZIR,I.A.F.,AUGUST 1987 
C16] a 
C17] a X AND HeTHE ARGUMENTS TO BE CONVOLVED (VECTORS) 
(181 A 
X 
تو‎ 1 
H 
Dua 1 
Yex LCV H 
Y 
9 12 104 1 
LONG 
T2 3 


AT ου. ου ο 15 14 15 16 17 18 19 20 21 22 23 24 25 26 27 2829 
30:31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 5253 
54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 7677 
78 79 80 81 82 83 84 85 86 87 88 89 90 


YeX LCV LONG 


Y 


3 8 14 20 26 32 38 44 50 56 62 68 


146 152 158 164 170 176 182 
254 260 266 272 278 284 290 
362 368 374 380 386 392 398 
470 476 482 488 494 500 506 


74 80 86 92 
188 194 200 
296 302 308 
404 410 416 
512 518 524 


95 


98 104 110 116 122 128 134 140 
206 212 218 224 230 236 242 248 
314 320 326 332 338 344 350 356 
422 428 434 440 446 452 458 464 
530 536 269 90 


(01 


۶701 7 
YeXl FLCV BALSOTO Ni: N2 Ne تر‎ 850 5 


(11 A 
(21 aA FAST 1-D LINEAR CONVOLUTION 
[31 A 
C4] 00+0 
[51 NlepX1 
(6] N2¢oH1 
C7] NeN1+N2-1 
C8] χΧΖεντχι 
C9] HZeNTH1 
[101 YeXZ FCC H2 
[111 >O 
[121 ñ 
C13] aA FAST LINEAR CONVOLUTION USING FCC FUNCTION 
C14] a 
[15] R8 Y.KATZIR,I.A.F.,AUGUST 1987 
C161 ñ 
C17] A X AND H-eTHE ARGUMENTS TO BE CONVOLVED (VECTORS) 
(181 A 
x 
35.3 
H 
3 2 
Yex FLCV H 
Y 
9 127102452 


YeX FLCV LONG 


PLEASE USE CC OR LCV FUNCTIONS. 
VALUE ERROR 
FLCVUC101 νο FCC ασ 


A 


96 


[O] 
EL] 
[2] 


[36] 
[371 
[38] 
[39] 
C40] 
C41] 
C42] 
C43] 
C44] 
C45] 


z 


VFFTIDJIV 
errr X:Kin:r:O0I0 

DIoeo 

2(re(2€6o0o0xX)^2€1T0x)18 
2(1€o0k)77-re25»kel0.5xox 

‘NOT A COMPLEX VECTOR’ 

+0 

>8 ,oxex, [70.510 

xe*Ck,2)0ox*2 

2(20zke 712 4 8 16 32 64 128 256 512 1024 2048 4096 8192 16384 3276 
>*Cllozex)l O O 

‘NOT A POWER-OF-2 ORDER’ 

--DEX 'z' 

Ζς(1 2 ,Ke2) 0x 

*kL 0 28 3 

>(2€0NC 'wt'2116 

2(k6ooxewt)T720 
ne0,-o6141xe20(0242n)xvin*4 
wtexe(iXo2)o(x,n2),n,-X 

220 


ne-/t21]z 

Ke€oze*/t21z 
2(4«k*€ooze€2,L[01]C-/L[1]nxkox2,[0.5]*/t1]nxKkoex)T7T19 
ne-/t21z 

ze*€/t21z 
nf;0511€-n[;0;1] 

nur 911]enEs:1] 

21-2 ,7 
ze9C-/2),Llneol]-/z 

LI 

ne “100z 

REIL: Ie=n[1:] 

xen-2 

n© ltoezez+n 
Κςο,ΦΙ19Υς2Ο(ΟΞΠ)ΧΙΠΞ2 
πο οσ1ροε,ςσ,-ε 
XeC(*zZxxr),[70.5]-7rxex 
Ze(Zz*x),2-xX 

30 


APL-STYLE FAST FOURIER TRANSFORM 
XeCAN BE EITHER REAL (A VECTOR) OR COMPLEX (2-ROW MATRIX) 
OPTIMIZED FOR SPEED 
Zz: THE OUTPUT IS COMPLEX (REAL AND IMAGE PARTS). 
SOURCE: PAUL PENFIELD,JR. 6/13/79 


ET 


[01 
[1] 
(21 
[3] 
C4] 


VIET 


See IEET ۶ ۸ 


DDDDD DVD 


Ώτοςο 
2(re(2€epox)^2€1T0ox)18 
2(l1€o0k)T7-re25»ket0.5xox 
'NOT A COMPLEX VECTOR' 
+0 

98,o0xe€ex,L 0.510 
xe9(k,2)px+4 


2(20xzke 712 4 8 16 32 64 128 256 512 1024 2048 4096 8192 16384 3276 


*(Cllozex)L 0 O 

"NOT A POWER-OF-2 ORDER’ 
~-OEX 'z' 

zeCl 2 ,ke2)pxen 

ək} O 28 23 

>(2€0NC 'wt')116 
>(Kkeppxewt)3T20 
nt0,-01l1xe20(022=n)x1in+4 
wte<exe(Ko2)o(x,n),n,-x 
220 
xesecko'xt::::::::5522 1:55 5G DONI 
ne-/t21z 

Κερζςν/{212 
$+(4<kepozez,10J(+/[1Inxko9x),10.5]-/[(1]Inxkp99-x3)T19 
ne-/t21z 

ze+/121z 

n. ;1:1l1e=n[:1l;11 

5 ۰۰0-1118 [7 
2e2,10J]n 
ze9(*/2),lne01-/z 

Vr 

ne {100z 

n[1;1€-nt1;] 

xez-n 

ne “ltezeze+n 
ke0,ó1lire20(0:+n)xin+2 
re(oz)ok,r,-r 
xe€(-—-Z#xxr),( 7”0.51+Zxxər 
Ze€(Z2*x),2-x 

0 


APL-STYLE INVERSE FOURIER TRANSFORM 

XeCAN BE EITHER REAL (A VECTOR) OR COMPLEX (2-ROW MATRIX) 
OPTIMIZED FOR SPEED 

SOURCE: PAUL PENFIELD,JR. 6/13/79 
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X 
۲1 30 IT11074l1l1131 


IED A 


Oo 


Oo 


Oo 


O O 


oo 


oo 


oo 


Oo 


oo 


oo 


oo 


w:o 


τες Y‏ ہد 


eo 
edo 
eo 
eo 
eo 
eo 
edo 
eio 
eo 
eio 
eio 
eo 
eio 
eo 
eio 
eno 


P 


XXexxLo 


BEEN 1111111111111 


29 


CO] 
EI] 
[2] 
[31 
C4] 
[5] 
[61 
CZI 
C8] 
91 
C10] 
(111 
[121 
[131 
C14] 
r 15,] 
[16] 
[17] 
(181 
[19] 
(201 


[ο] 


(101 


VPSECOIV 
INWeZP PSE DATA;UIO;N;XN;XW 

A 

A COMPUTING THE POWER SPECTRUM ESTIMATION (PERIODOGRAM) 
A 

DIO<O 

Ne(pDATA)+ZP 

Ne2*f (20N) 

XNeDATA, ((N-(pDATA))p0) 

XWeFFT XN 

INWe (+4 (XW*2) )=N 

30 


COMPUTING THE POWER SPECTRUM ESTIMATION (PERIODOGRAM) USING FFT 
FUNCTION 


Y.KATZIR ,I.A.F. AUGUST 1987 
DATAFTHE DATA VECTOR (1-D) 


ZPeLENGTH OF ZERO PADDING FOR SMOOTH PLOT 
INW: THE OUTPUT-SX(0) 


DDDDDDDDODOD BD 


VPSEBLOUIV 

INBWeSEG PSEB DATA; DIO: M:¡MM:XW;XN:L:;X 
A 

A COMPUTING THE POWER SPECTRUM ESTIMATION USING BARTLETT WINDOW 
A 

OI0<0 

Mel ( (eDATA)=SEG) 

XNe(M,SEG)p0 

ZP+*256-SEG 

MMe(1lo0XN3+2ZP 

XWe(M,MM)p0 

LeKeo 

LOOP1: XN(CK3; LSEGJ]©DATAL (SEGxK) + iSEG] 

$ (M>KeK+1)/LOOP1 
LOOP2:XW[L;1€ZP PSE XN(L;] 

> (M>LeL+1)/LOOP2 

INBWe(+2XW)=+M 

30 


Y. KATZTIR,I.A.F. AUGUST 1987 


SEGeLENGTH OF SEGMENTS 

DATA+THE DATA VECTOR (I-D) 

USING IHE PSE FUNCTION (+FFT). : 

THE FUNCTION PADS THE SEGMENT WITH ZEROS UP TO 256 POINTS. 


DDDDDDDD 


100 


vs Y 


e 
wá 


Pilot oF 


v» u) u>) 
u) u u) 
u? Ud) η 
u) U? u»5 
u? u? u? 
u? u? u? 
Www 


150 


u) u) u? 
u) U) U) 
u? u? WwW 
u) u? u? 
v) uv? u 
u) u) u) 
wm Uy 4) 
ο uv» v) 
vo wm 
uy UD v) 
v) Ww» v» 
mw س١‎ 
uU» v» v») 
u) U) u) 
u) u u? 
v) wv» y) 


= Y 


== 
9 مک 


u) u) u) u) 
u? Vv) UD Uu» 


Plot of 


Www Uu) 


101 


ud 
u) 
u? 


ھا 
u‏ 
uy‏ 


u? 
u 


2 





Ye240 PSE X 
€128¢ 
Yev 128 


e 
*- 


uU» v) u?) u> 
mm u) u?) u) 
u) U) u?) یا‎ 
V» v) UD Uu) 
u) U) Ww uy) 
min uv» Ww 
wu» UN UN wm 
u) u) UY) v) 
ww Ww 
u> U) u) u? 
v» YU) YN m 


( 4 segment of 32 points. 2 





Ze32 PSEB X 


£25 $ 


eZ 
2707 


uv» 


z N 


u> N 


| 7۳ ؟ 
XZeXSQRT XR‏ 01[ 


[1] A 
[2] A COMPLEX SQUARE ROOT (XR*0.5) 
CIEGA 


C4] X2€C0.5-0^.»XROxC2xCC*7XRXxXR2*0.5)*1 10 *.xXR)*0.5 
(5] XZe(0A.> 1 0 7XR0eX2,[010-0.51CO0 1 *.xXR)*2xX2-*X270 


X 
0.5 0.866 


XSQRT X 
0.86602 0.49999 


VXTIMEL(UIV 
CO] XZeXL XTIME XR 
(11 A 
[21 a COMPLEX TIMES (XLXXR) 
(31 a 


C4] XZe((1 0 7XLOx 1 O 7XR2-CO 1 7XLOx O 1 7XR 
(51 XZex2,L0I01CC1 O 7XLOx O 1 7XRO*CO 1 7XLOx 1 O 7XR 


X 
0.5 0.866 


X XTIME X 
“0.49996 0.1866 
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۳۸۹7 
[01 ZeXMAGN XR 


[11 a 
(21 a COMPLEX MAGNITUDE (!IXR) --RETURN IS REAL 
(31 A 


C4] Z€(*7XR*2)*0.5 


X 
0.5 0.866 
XMAGN X 
0.99998 
VXPHAS(OIV 
[01 ZeXPHAS XR 
[11 a 
(21 a COMPLEX PHASE ANGLE --RETURN IS REAL 
[31 A 


C4] Ze O 1 +.XXR 
[5] Ze(xZ2+2=0)x(0v.f*XRIxX"20(1 O +.xXRI>(0A.=XR)+(+2%*XRXXR)*0.5 


x 
0.5 0.866 


XPHAS X 
1.0472 
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VXCONJCOIV 
CO] XZeXCONJ XR 


(11 Rñ 
(21 A COMPLEX CONJUGATE (+XR) 
[3] A 


C4] X2ecl 0 τικ) πο O 1 #XR 


X 
0.5 0.866 
XCONJ X 
0.5 70.866 
ΘΧΕΧΡΟΓΩ͂ΙΥ 
[0] XZeXEXPO XR 
[11 ‘a 
(21 a COMPLEX EXPONENTIAL (*XR) 
[3] aA 


C4] XZeCCoXR)o* 10 XXRIX 2 1 9.0 O 1 *.xXR 


X 
0.5 0.866 


XEXPO X 
1.0682 1.2559 
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VELCOIV 


(0]  Yepr 
[L1] Yeol 
pI 
3.1416 
"900 ۳ 
CO] Y*cos x 
CT] Y#20X 
COS (PI+4) 
0.70711 
VSINLDJIV 
CO] YeSIN x 
[1] Yelox 
SIN (PI+2) 
1 
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JLOAD O EC4440 
O EC4440 SAVED 2/23/1988 14:26:16 


| NAVAL POSTGRADUATE SCHOOL MONTEREY, CALIFORNIA | 


SP d ٭-‎ eS “ms “sm “s A A amm A E A “mm “=s کے ہے‎ A A Qz e تو چ‎ m c m “GE “sm A ED AE ےے>‎ “m “m “m < “m که‎ A A س س‎ me ۱ 


ECE-4440 WORKSPACE 


۱ ۱ 

I ۱ 

۱ THESIS BY ۱ 

۱ ¥.XATZIR.I.A.F. ۱ 

۱ ADVISOR: PROF. C.W. THERRIEN. » | 

| VERSION 1.0 SEPTEMBER 1987 | 

JFNS 

BARTLETTC BARTLETTR cc CC2D CHEB COEFF Cos FCC 
FCC2D EFT EFTI2D FLCV HAMMINGC HAMMINGR HANNINGC 
HANNINGR LEFT IFFT2D  ILPFILT INFORMATION LCV MCCLEL MEAN 
PE PROTFILT RECTANR SHIFT SHIFT2D SIN TRANSFNC UNSHIFT 
UNSHIFT2D XCONJ XEXPO XMAGN XPHAS XSQRT XTIME 
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(01 
C11 
[2] 
C3] 
(41 
E51 
(61 
(71 
C8] 
(91 
(10] 
T11] 
[121 
1 
[141 
[151 
C16] 
t17] 
[18] 
51 
[201 
[21] 
022] 
C23] 


۹67 
οι υπ τι το 000:072 2 72 
A 

A 2-D CIRCULAR CONVOLUTION 
A 

>((A/(oeoX)=oHB)A2=Q=poX) /BEGIN 
De’ IMPROPER ARGUMENTS’ 

20 

22601:0700 

Nleltox 

N2¢110X% 

Y*(N2,N1)p0 

Í1+*0 
LOOP1:T1€C-CI1*1))99H 

1250 
LOOP2:T2€(-(I2*1))0eT1 
Y(12;11]e+/+/XxT2 
-(N2»I2€12*1)/LOOP2 
»-(N1»I1€I1*1)/LOOP1l 

»0 


THE INPUTS NEED TO BE 2-D AND THE SAME SIZE. 


Y.KATZIR 1.A.F.,SEPTEMBER 1987 


DDDD D 
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2 2 10 
1100 
0000 
0000 
H 
1100 
1100 
1000 
0000 
Yex CC2D H 
Y 
2431 
3641 
3420 
1100 
21 


l2 3 45 67.89 10 11 12 13 14.715 I6 


X CC2D 21 
IMPROPER ARGUMENTS 


22 


مم س 
eN‏ 


X CC2D 22 
IMPROPER ARGUMENTS 
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[01 

El] 

[21 

31 

C4] 

[51 

(61 

[7] 

(81 

(91 

C10] 
E] 
[121 
EX3] 
(14] 
[15] 
C16] 
[171 
(181 
(19] 
(20] 
C21] 
(221 
[23] 
(241 
(251 
C26] 
p27] 
C28] 
E29] 
(301 
ΠΠ 
C32] 
G33] 
C34] 
55: 


VFFT2DLO1V 
YeFFT2D X;0IO;NROW;NCOL;1;J;G 


A 
A 2-DIMENSIONAL FFT 
A 
OL0+0 
A CONVERT REAL ARRAY TO COMPLEX ARRAY 
I#p (eX) 
+9(I=3) / CONTIN 
NROW#p XC 31] 
NCOLepXC13 4 
Ge(2,NROW,NCOL) p0 
GL0;;lex 
XeG 
A COMPUTATION OF FFT 
CONTIN: NROWep XC 1331] 
NCOL*0oX[1;1;1] 
Gexxo 
Yes 
IeJ€0 
COLLOOP:G[;;Jl]€FFT X([;;J1 
»(NCOL>JeJ+1)/COLLOOP 
ROWLOOP:Y[;I;]eFFT GL[;I;:] 
3 (NROW>TI*eI+1) /ROWLOOP 
A IF OUTPUT ARRAY IS REAL VALUE,2ND PLANE IS DELETED. 
το ον στις] 
ο ο 1ο έτος] 
I«€9/5/*5/CG£*Y) 
~(I=0)/STPNOW 
20 >` 
STPNOW: Ye C1,NROW,NCOL) TY 
Y#CNROW,NCOL) oY 


FFT FUNCTION HAVE BEEN USED. 


Y.KATZIR,I.A.P. SEPTEMBER 1987 
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C0) 
EI 
[2j 
[3] 
(41 
(541 
C6] 
[71 
C8] 


C10] 
[IlI 
[121 
[131 
[141 
[15] 
C16] 
[17] 
C18] 
[19] 
C20] 
[211 
[221 
C23] 
(241 
[251 
C26] 
C271 
C28] 
C29] 
C30] 
[311 
C32] 
[331 
[341 
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SHIFTING 2-D PERIODIC DATA TO AN INTERVAL CENTERED AT THE ORIGIN 
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A IDEAL LP FILTER (CIRCULAR BASE) FREQUENCY RESPONSE 
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We(P,PJ)p0 

IeJ+0 

PBe(P=2)-1 

RBe(PxXR)+(2xPTI) 
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ZERO: 9( P>Je¢J+1)/LOOP 
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ReTHE RADIUS OF THE BASE (GAIN=1) 
PeTHE DIMENSION OF THE MATRIX 
W: THE IDEAL FREQUENCY RESPONSE 
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II+WINDOW BASE DIMENSION 
NeMATRIX DIMENSION 
WR: THE WINDOW MATRIX 
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IIeWINDOW BASE DIMENSION 
NeMATRIX DIMENSION 
WR: THE WINDOW MATRIX 


Y.KATZIR,I.A.P.,SEPTEMBER 1987 


ο ο ο ο ο ο Ὁ 


123 





2) 


az 


r ba 


HANNING WINDOW ‘Peetangyls 


10 


جع وو ا سو ا 


URe[T HANNINGR N 
9 ΠΡΟΚΣ 115 


[ieli 
He32 


"-*9o00*29:9 9920*9.299*20080*090!752*9 


124 


ο ο λα ΣΕ ٥ 
""*9*200202:9009*09.08922*9898*:80.09*224€4 


سی سس ری رٹ و رز رز و روز رڈ رز زور رر .دی 


SPcSceerereesereatsersensesat jars 


SSCeSrcaaterettiasF@otes®Sreases co 





νο. Pos ΚΝ. [ISZITENIT 


eerer'e ereten 0. دمح بوڈ درم‎ eens 
9* 898180 @b 000... 


contour plot 
Be 


ve‏ وو او وا یہہ 


eetec@eerces®e — —  @8@2e¢¢n08 
HO... e vere... erecto pa. 
SeeSe Serer eetarsaSeSsOereSsStiarss 
SSeS SSese SF eSteSeat?Ge Ge GeO sare 


PerSeeSe SeF St یی رر‎ - ٔ ٔ  -- 


ο ΕΕ ΕΣ Ἡ 


SetrStraGeasFvsesF®astetSetasar seen 


VHANNINGC(OIV 

WCeII HANNINGC N;DIO;I;J;PR;RB;RR 

a 

R CIRCULAR BASE HANNING WINDOW 

a 
010+0 
WC+*(N,NJpO0 
2 0ے نے‎ 
PRe(N+2)-1 
RBe(IIxN)=(2xPI) 

LOOP: RRe(((I-PRI*2)+((J-PRI*2)I*0.5 
+» (RB<RR)/CONT i 
WC[J;¡IJ“0.5x(1+(COSC(CPIXRRI+>RB))) 

CONT: 3(N>J*eJ+1)/LOOP 
J+0 
»*(N>I¢I+1)/LOOP 
20 


II+THE RADIUS OF THE BASE 
NeTHE DIMENSION OF THE MATRIX 
WC: THE HANNING WINDOW 


DDDDDDD 


Y.KATZIR,I.A.F., OCTOBER 1987 


125 


BLONET SJ! 


N€ 32 


UCeI] HANHIHGC H 


0 ΠΡΟΧΕ 116 





contour plot 


L 
9 


رر رْ یر شر پٹھڑيے. 99*:9:09»09590:9.005.9140 
Φ»ο» 89»‏ 9599551555 و و رز رر یرب 
AI 000 "9*8‏ 


"9*:9:0906:025:221 ***9*o:220c2c2: ... 
»5.:5.»9945» رر ہب‎ َ ىىى۷ٰٰىئ۷۷س‪‪م-٘<0٥۵‎ 
dor ..ro.s *"9*9*9os9ve00s0c050ss2^4 
?9*9s:909:0c00520y e (eis ioco 

LEE E ET SO Stevresseta‏ کر رر ور رب 
ےمج روہ وھ .%6010.090.:0:9200 
E 6409000015‏ 9919:90:»9»90 

. 


. 9 
LANE . o: ee tt... * 
.. . - s > = e e n.. 
9919199659: 999 , ort Hr: 106 0.. 
a φ 
@ ee ..... 9: 9**9 docto 8 υσδουσ.ο 
5916919929» bono o 5999 959521950 


..... 9595594 Sere t... Pro. .... 


. . m 
...... eso: 90606 5 ος Pret tore 
90’ 

5.505 05509505 αυ ν 
. 





126 


(0] 
CI] 
(2] 
US] 


[221 
ΓΞ J 
C24] 


VBARTLETTR(O1V 
WReII BARTLETTR N;U0IO;HW1;HW2;J;I;ROW;COL;S 
A 
A RECTANGULAR BASE BARTLETT WINDOW 
A 

OI0<0 

WRe(N,NJ 00 

EW1+*HW2*IIp0 

Se(II-1)+2 

ROW+*0 
ROWLOOP: COL+0 

I#ROW+(2xS) 
HWl1[ROW1e€1-ClCROW-S))-4S 
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APPENDIX D 
COMPUTER ASSIGNMENTS AND SOLUTIONS 


This appendix provides a few samples of computer assignments and related 
solutions using the software package. These are in addition to the examples in 
Chapter III. The computer assignments in EC 3400 were taken from FIRST PRIN- 
CIPLES OF DISCRETE SYSTEMS AND DIGITAL SIGNAL PROCESSING by 
Professors R. D. Strum and D. E. Kirk [4]. The computer assignments in EC 3410 


and EC 4440 were taken from Professor C. W. Therrien’s homework assignments. 


EC 3400 COMPUTER ASSIGNMENT [4] 


The purpose of a bandstop filter is to remove the effects of the vibrations in 
the frequency range 6 Hz to 15 Hz from a control loop. The sampling frequency is 


fs=100 Hz. 


a. Using the Fourier Series design procedure, calculate the filter coefficients with 
a rectangular window. 


. Plot the magnitude portion of the frequency response for a filter having thirty- 
one coefficients. 





. Repeat part (b) using a von-Hann window. 


SOLUTION 
1. By using the function DIGFREQ the digital upper and lower frequencies are 


calculated (Figure D.1). Figure D.2 describes the ideal frequency response to 
be generated. 
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Figure D.1. Calculating the Digital Frequency Using the Function DIGFREQ. 
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Figure D.2. The Ideal Frequency Response. 


2. By using the function BSCOEFF 31 causal impulse response coefficients are 
calculated (See Figures D.3 and D.4). 

3. The function HANNING generates 31 von-Hann window samples that are 
multiplied by the filter coefficients. The result is windowed coefficients (See 
Figures D.5 and D.6). 

4. By using the function FREQRES the frequency responses with and without 
window are generated. The following shows the APL commands which gener- 
ate the frequency responses, and Figures D.7 and D.8 show the plotted results. 


HBS—128 FREQRES hBS 
HBSW—128 FREQRES hBSW 
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hBSe31 BSCOEFF TUL 


ABS 
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0.82 70.14034 70.042416 0.063217 0.12619 0.12421 0.07206 7.8548E' 3 
70.032854 70.037409 70.01871 71.0218E^3 78.2842E^4 “0.016485 70.03 
70.033694 


Figure D.3. Calculating the Impulse Response Coefficients. 
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Figure D.4. The Impulse Response Coefficients. 
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Figure D.5. Generating the Hanning Window. 





Figure D.6. The Hanning Window. 
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Figure D.7. The Frequency Response Magnitude. 
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Figure D.8. The Frequency Response Magnitude with Window. 
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EC 3400 COMPUTER ASSIGNMENT [4] 
A nonrecursive bandpass filter is to be designed using the IDFT approach. 


The ideal characteristic is given below: 


5F 


— TT -m = -- - 


6 2 


Use 255 frequency samples and find 51 filter coefficients using: 
a. A rectangular window 


b. A Hamming window 


c. Plot the frequency response for the filters of parts (a) and (b). 





SOLUTION 


1. By using the function IDEALF and choosing 256 samples, the ideal bandpass 
frequency response is generated (See Figure D.9). 


2. By using the function FCOEFF 51 coefficients are calculated (See Figure 
D.10). 


3. The function HAMMING generates 51 Hamming window samples that are 
multiplied by the filter coefficients. The result is the windowed coefficients. 
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Figure D.9. Generating the Ideal Frequency Response Using 255 Sam- 
ples. 


he51 PCOEFF H 


h 
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0.036326 70.069496 0.025239 0.024298 5.1718E^3 70.029088 ^0.077918 
0.21179 70.12868 70.16131 0.32813 70.16131 “0.12868 0.21179 “0.077 
70.029088 S.1718E”3 0.024298 0.025239 70.069496 0.036326 0.011158 
75.0631bE^3 70.013214 70.010498 0.040397 70.025658 “S.S5166E”3 4.885E”3 
8.7724E^3 4.7036E^3 70.027498 0.020547 2.582E^3 74.642E^3 “6.2133E73 


Figure D.10. Calculating the Impulse Response Coefficients. 
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Figure D.11. Generating the Hamming Window and Calculating the 
Frequency Response. 
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Figure D.12. The Frequency Response Plots. 
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Figure D.13. The Frequency Response in dB 
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Figure D.14. The Frequency Response Phase. 
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EC 3410 COMPUTER ASSIGNMENT [5] 






On a floppy disk you will find 4 data sets called 









SLO.DAT 


SL1.DAT 


SL2.DAT 






SL3.DAT 






Each of these represents a random signal with 512 time points. The format of each 






file is identical. You can look at the file with an editor. The first line contains the 


number 512 displaying the number of points in the data set. On the remaining 







lines, having four numbers per line, you will find the floating point values of the 


signal. 






1. Transfer these files to your APL WS using the UTILITY function GETDATA. 





2. Plot each of the signals versus time and submit the plots with your assign- 
ment. Tell whether you think each of these signals seems to be more-or-less 
uncorrelated, positively correlated, or negatively correlated. 









. Compute the mean of each signal by averaging in time and write down your 
result for each signal. It is reasonable based on your plots? 





. Subtract the mean from each signal and compute and plot the sample corre- 
lation function. Attach the plots to your assignment. Was your guess about 
the correlation of these signals in part 2 correct? 





SOLUTION 


1. Using the function GETDATA and the following statements, the data sets 
SLO, SL1, and SL3 are transferred to the APL workspace: 


SLO<—GETDATA ’SLO.DAT’ 
SL1c—GETDATA ’SL1.DAT 
SL2—GETDATA 'SL2.DAT 
SL3—GETDATA *SL3.DAT 
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2. By using the X-Y plotting procedure in STATGRAPHICS the data plots are 
generated. (See Figures D.15, D.16, D.17, and D.18.) Another kind of plot 
can be obtained by using the STATGRAPHICS Time Series Procedure. The 
results are shown in Figures D.19, D.20, D.21, and D.22. 
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Figure D.16. SL1 Data Set. 
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Figure D.17. SL2 Data Set. 


Urn 





Figure D.18. SL3 Data Set. 


3. From the plotting of each data set, it can be seen that: 


a. SLO changes randomly with little consistency, this is characteristic of a 
white noise sequence. 


b. SLI varies more slowly and stays stable for a while with positive values 
and then changes more-or-less smoothly to negative values. This indicates 
that it is a positively correlated data set. 


c. SL2 is very nearly like a white noise sequence. 
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Figure D.19. SLU Data Set Using Time Series Procedure. 
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Figure D.20. SL1 Data Set Using Time Series Procedure. 


d. SL3 changes at almost every step from negative to positive. This is char- 
acteristic of a negatively correlated data set. 


4. Calculate the mean of each data set using the MEAN function. Then subtract 
the mean of each data set from the data set in order to achieve a zero mean data 
set. Using the following instructions, calculate the first 150 autocorrelation 
lags for each data set: 


Rxx1«—150 SAC? SLO 
Rxx2—150 SACF SL1 
Rxx3<—150 SACF SL2 
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Figure D.22. SL3 Data Set Using Time Series Procedure. 


Rxx4—150 SACF SL3 


Using the X-Y plotting procedure in STATGRAPHICS the estimated auto- 
correlation is plotted as shown in Figures D.23 to D.26. 

Another method for producing the autocorrelation values and plot them 1s to 
use the autocorrelation procedure in the TIME SERIES ANALYSIS section from 
the STATGRAPHICS Menu. The results are shown in Figures D.27 to D.30. 
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Figure D.24. Sample Autocorrelation Function of Positive Correlated 


Set. 
Looking at the plots verifies that our characterization of the data sets in para- 


graph 3 was correct. 
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Figure D.26. Sample Autocorrelation Function of Negative Correlated 


Data Set. 
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Figure D.28. Estimated Autocorrelation Using Time Series Procedure. 


158 


(«2389 707900 





o 4o So 120 180 
T 


| | 1 ! H! null] 


00 


οσο” O 0) 





O LO So 420 160 
las 


Figure D.30. Estimated Autocorrelation Using Time Series Procedure. 
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EC 3410 COMPUTER ASSIGNMENT [5] 

Compute power spectrum estimates for the data sets SLO, SL1, SL2, and SL3 
using a rectangular window. In all cases zero-pad your data appropriately before 
giving it to the FFT routine so that you will have computed enough points of the 
spectrum to give a smooth plot. (Plotting say 256 points of the spectrum should 
be sufficient). Assume the data has been sampled at a spacing of T = 0.1 ms. 
Label the frequency axis in hertz. Since the signals are real, it is only necessary to 


plot the positive frequencies. 
a. Use only one 16-point segment (the first 16 points) of each data set. 
b. Use 32 16-point segments. 
c. Use 16 32-point segments. 
d. Use 4 128-point segments. 


What do you make of the results? 


SOLUTION 


1. The first step is to extract the first 16 samples of each data set and compute 
the power spectrum estimate using the function PSE. In order to obtain a 
smooth plot the data sets are zero padded with 240 points. (See Figure D.31.) 
The power spectrum estimate will have 250 points as shown on the plots of 
Figures D.32 to D.35. 


2. By using the function PSEB (power spectrum estimation with Bartlett win- 
dow), we estimate the periodogram of the SLO data set using 30 segments of 
17 points, 16 segments of 32 points, and 4 segments of 128 points as seen in 
Figure D.36. Figure D.37 is the periodogram plots for each case. The above 
process is repeated with the data sets SL1, SL2, and SL3. See Figures D.38 
to D.43. 


3. From these results we conclude that the first data set has a spectrum like 
white noise. The second and third spectra are lowpass with bandwidths of 
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Figure D.31. Calculating the Power Spectrum Estimates Using the 
Function PSE. 


approximately 1 kHz and 1.2 kHz. The fourth spectrum is highpass with 
bandwidth of approximately 1 Khz Hz. It seems that the best estimate is 
obtained by using 16 segments of 32 points each. 
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Figure D.32. The PSE of SLO Data Set Using 16 Points. 
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Figure D.33. The PSE of SL1 Data Set Using 16 Points. 
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Figure D.34. The PSE of SL2 Data Set Using 16 Points. 
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Figure D.35. The PSE of SL3 Data Set Using 16 Points. 
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Figure D.36. Calculating Periodograms of SLO Using Bartlett Win- 
dows. 
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Figure D.37. Periodograms of Different Windows. 
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Figure D.38. Calculating Periodograms of SL1 Using Bartlett Win- 
dows. 
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Figure D.39. Periodograms of Different Windows. 
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Figure D.40. Calculating Periodograms of SL2 Using Bartlett Win- 
dows. 
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Figure D.41. Periodograms of Different Windows. 


169 


5531517 PSEB SL3 


SL31+1281SL31 


SL32€32 PSEB SL3 


SL32©128TSL32 


55335128 Ρ5ΕΒ 513 


SL33©128TSL33 


Figure D.42. Calculating Periodograms of SL3 Using Bartlett Win- 
dows. 
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Figure D.43. Periodograms of Different Windows. 
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EC 4440 COMPUTER ASSIGNMENT-FILTER DESIGN [6] 

This problem deals with the design of a lowpass filter by the windowing 
method. Use APL to do this problem. Use STATGRAPHICS to generate 3-D 
and contour plots. 


1. Start with the ideal frequency response defined by the following equation: 


1 , wy +wz < (0.47)? 
I e 1 2 
(2سا,1س)‎ | O , otherwise —-T<w],w ST 


Construct an APL array with its values. Note that for compatibility with the 
FFT, the range on w, and wa should be zero to 2 rather than —r to +1. Use 
a 32 by 32 point array. | 


2. Find the ideal impulse response :(n1,n2) and plot it. 


3. Choose one of the 1-D windows and plot it. Define a 2-D rectangular window 
from it. Multiply 7(nl1,n2) by the window to get h(n1,n2). The window 
should by 11 by 11 points in size. 


4. Use the Fourier Transform (FFT) of h(n1,n2) to get the filter frequency re- 
sponse. Plot this using both 3-D plotting and contour plotting. 


SOLUTION 


1. By using the function ILPFILT (Ideal LP filter) samples of the frequency 
response of an ideal lowpass filter is generated (32 by 32 point array). The 
plots achieved by using the surface plotting procedure in STATGRAPHICS 


can be seen in Figure D.44. 


2. By using the IFFT2D function (Inverse 2-D Fast Fourier Transform) the ideal 
impulse response is founded and plotted in Figure D.45. 


3. By using the function HAMMINGR we generate a 2-D HAMMING window 
with rectangular base. See Figure D.46. 


4. Figure D.47 shows the product of the HAMMING window and the ideal im- 


pulse response. 
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9. By using the function FFT2D (2-D Fast Fourier lransform) the frequency 
response of the windowed filter is found, and by using the function XMAGN 
the magnitude of the result is calculated and plotted as in Figures D.48 and 
D.49 (contour plot). Note the use of the function SHIFT2D to show the 
frequency response in the interval 0 to 27. 
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Figure D.44. Calculating and Plotting the Ideal Frequency Response. 
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Figure D.45. Calculating and Plotting the Ideal Impulse Response. 
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Figure D.46. Generating and Plotting the Hamming Window. 
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Figure D.47. The Product of the Hamming Window and the Ideal Im- 
pulse Response. 
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Figure D.48. The Design Frequency Response. 
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Figure D.49. The Shifted Design Frequency Response. 
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